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High harmonic generation by relativistically intense laser pulses from overdense plasma layers 
is surveyed. High harmonics are generated in form of (sub-)attosecond pulses when the plasma 
surface rebounds towards the observer with relativistic velocity. Different cases are considered. The 
"relativistically oscillating mirror" (ROM) model, describing the most typical case, is analyzed in 
detail. The resulting harmonic spectrum is usually a power law with the exponent -8/3 jl], but 
possible exceptions due to "higher order 7-spikes" are considered. It is shown that under certain 
conditions, ultra-dense electron nanobunches can be formed at plasma surface that emit coherent 
synchrotron radiation. The resulting spectrum is much flatter and leads to the formation of a 
giant attosecond pulse in the reflected radiation. The harmonics radiation is also considered in time 
domain, where they form a train of attosecond pulses. It is characterized and a possibility to select a 
single attosecond pulse via polarization gating is described. Further, the line structure in relativistic 
^ harmonic spectra is analyzed. It is shown that the harmonics have an intrinsic chirp and it can be 

^ O responsible for experimentally observed spectral modulations. Finally, high harmonic generation is 

considered in realistic three-dimensional geometry. It is shown that free space diffraction can act as 
a high pass filter, altering the spectrum and the reflected fleld structure. The high harmonics tend 
to be self-focused by the reflecting surface. This leads to a natural angular divergence as well as to 
field boost at the focal position. Coherently focusing the harmonics using an optimized geometry 
may result in a significantly higher field than the field of the driving laser. 
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I. INTRODUCTION 

in 

Qh The tremendous progress in the femtosecond laser technology leads us to the question, if it is possible to create 
even shorter pulses, with durations in the attosecond or, may be, zeptosecond time range. 

O One, already well-established way to produce light pulses of attosecond scale duration is the generation of high 

c/3 order harmonics by the non-linear interaction of gas atoms with lasers of intensities close to the ionization threshold 
[2]. This mechanism is however limited to not-too-high laser intensities and a relatively low efficiency. Production of 
coherent attosecond pulses of higher energy could open up the way to entirely new methods of attosecond research 
I such as XUV pump-probe spectroscopy [31 13]. The currently most promising way towards more intense attosecond 

pulses is the generation of high order harmonics (HHG) at solid density plasma surfaces. 
^-H For a complete understanding of this attosecond pulse generation scheme, it is necessary to study three stages in 

^ the laser-plasma interaction process: 

m 

Plasma Formation. Before the main laser pulse hits the solid target surface, the pedestal of the pulse already ionizes 
' . ^ it and turns it into a plasma. The plasma then thermally expands and at the same time is pushed inside by the 

laser ponderomotive potential. Depending on the contrast ratio of the laser system and the exact structure of 
the pre-pulse, this may yield very different surface density profiles. These processes are well understood today 
and can reliably be simulated by hydro-codes such as Multi-FS [5 . 
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Harmonics Generation. The second stage is the harmonics generation itself. It happens during the interaction of 
the main laser pulse with the pre-formed plasma density gradient. If the laser pulse duration is in the order 
of just a few ten femtoseconds or below, the motion of the ions during this period can be neglected, and the 
interaction takes place between the laser electromagnetic fields and the plasma electrons. 

Diffraction. After the radiation has been emitted from the surface, it will propagate through space. Due to the 
extremely broad spectrum of the emitted radiation and its coherent phase properties, it is well worth to take a 
closer look at its diffraction and focusing behaviour. 

The paper is organized as follows. First of all, in section |TT] we provide a short historical overview of high harmonic 
generation from overdense plasmas in relativistic regime. 
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In section [nH we examine the theory of relativistic HHG. Different models are discussed and spectra are analytically 
derived from the models via asymptotic analysis. This method has the merit of yielding "universal" spectra in the 
sense that these are independent on details of the electron motion, but only on some basic properties around the so- 
called 7-spikes. The most typical case is the power law decay I oc a;"*/'^ [B], runing up to a critical harmonic number 
ujr/tdo ^ 7'^, where 7 is the relativistic 7— factor of the oscillating plasma surface. At higher harmonic numbers, 
the spectrum rolls over into an exponential decay. Although the spectrum j6] holds in most cases, special and even 
more efficient regimes are possible. We found that the most efficient regime of single pulse HHG is governed by the 
formation of highly compressed electron nanobunches in front of the surface and results in a slowly decaying spectral 
power law with an exponent of 6/5 - instead of 8/3. The analytical results are substantiated by numerical simulations. 

In Section HY] we consider properties of this novel radiation source. Particularly, we discuss the structure of the 
attosecond pulses. 

Section |V] deals with the line structure of the high harmonics spectra. It is demonstrated that the generated 
"harmonic" lines can appear severely modulated due to the unequal spacing inside the attosecond pulse train. Thus, 
these spectral modulations do not imply incoherent noise, but rather are the result of a natural frequency chirp of 
the harmonics. We show that the line structure contains information about the motion of the surface plasma on the 
femtosecond timescale. 

Finally, section VI treats the surface HHG in a realistic 3D geometry. Here, diffraction takes on an important role 
shaping the radiation spectra. With carefully designed surfaces or laser pulses we may harness diffraction as a sort 
of spatial spectral filter. With a well designed focusing geometry, it should even be possible to focus the harmonics 
coherently in both space and time, yielding unprecedented intensities that exceed the intensity of the laser itself by 
more than a thousand times. 



II. HISTORICAL OVERVIEW 



The basic idea for HHG at overdense plasma surfaces has been around for almost thirty years now and endured 
several generations of high power lasers. In this work, we focus on the most efficient, highly relativistic regime. Before 
we move on to the actual study, it is worth to have a brief look at the history of the topic. 

The first observation of high harmonics from plasma surfaces was reported from the Los Alamos Scientific Laboratory 
[II H] in 1981. At that time, huge CO2 lasers were used at nanosecond pulse duration and the observed radiation 
was non-coherently emitted into the whole half-space in front of the target. A theoretical explanation for this first 
observation was given in Ref. The spectrum extended up to a sharp cutoff, which was found to be the plasma 
frequency corresponding to the maximum electron density. Therefore, non-linear collective plasma behaviour could be 
identified as the source of the harmonics. In the strongly inhomogeneous plasma, laser light was resonantly converted 
to plasma oscillations, which in turn produced harmonics by sum frequency mixing with the laser light. 

For some time then, it became silent around surface HHG, but interest rose again, when the CPA technique, 
invented in 1985 by Strickland and Mourou [lO,, revolutionized ultraintense laser science in the 1990's. With the 
newly possible fs-duration, multi-TW pulses, HHG entered an entirely new regime |111 112|. Because of the much 
shorter pulse duration, the plasma surface is not destroyed by the pulse and the harmonics are cleanly emitted around 
the specular direction along with the reflected fundamental [13]. 

In the mid-nineties, there were first theoretical reports about a novel HHG mechanism based on a non-linearity 
of purely relativistic origin, providing a source for harmonics without the limitation of a strict cutoff at the plasma 
shelf density |14H17| . The mechanism could roughly be described by a simple model, now commonly termed the 
"relativistically oscillating mirror" (ROM) [15, 16 . However, for the time being, lasers were still not strong enough 
to unambiguously demonstrate the relativistic effect in distinction to the non-relativistic plasma non-linearities. 

In the first decade of the new millennium, theory of surface HHG made further substantial advances. It was found, 
that for fs-laser systems the harmonics due to the plasma non-linearity were much stronger than could be expected 
from the old theory. This was attributed to so-called Brunei electrons [TS] that trigger the plasma oscillations instead 
of the evanescent laser field, leading to "coherent wake emission" (CWE) [19]. The ROM model was put on a solid 
basis by Baeva, Gordienko and Pukhov (BGP) who managed to calculate a universal spectral envelope (/ cx cj^*/"^) 
for the model by means of asymptotic analysis. The BGP theory takes fully into account the surface acceleration, 
leading to a smooth spectral cutoff at a frequency scaling as Wc c< 7^, comparable to synchrotron radiation, and not 
just proportional to 7^ like the Doppler frequency upshift at a constantly moving mirror. Around the same time, 
experiments were first able to unambiguously demonstrate the relativistic mechanism and confirmed the spectrum 
obtained in the refined ROM model, see Ref. [20J. 

In 2010, an der Briigge and Pukhov |2T] discovered another mechanism based on the relativistic non-linearity. They 
found out, that for certain combinations of parameters, extremely dense and narrow electron bunches may form at the 
surface. In this exciting regime, not even the basic boundary condition of the ROM model is valid and the frequency 
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upconversion process can be much more efficient than predicted by the model, 
coherent synchrotron emission (CSE) from the electron "nanobunches". 



The radiation is then described as 



III. GENERATION PROCESSES AND MODELS 



We discuss the theory of HHG at surface plasmas, with a focus on the highly relativistic regime ao ^ 1. 

In subsection |III A[ we start by summarizing the theoretical framework all models of the interaction are based on. 
Once having the equations written down, it is straightforward to derive some selection rules for the parity and the 
polarization of the generated harmonics. 

In subsection |III B[ the ROM model is discussed. We are going to see that this reputedly well-known model still 
bears some surprises, with respect to both its foundation and the spectra that can be derived from it. 

After this, we consider a variation of the ROM boundary condition that was suggested in Ref. [25] in conjunction 
with a two-pulse-scheme (subsection IIICl. We show, that this condition represents a "totally reflecting oscillating 
mirror" (TROM) with negligible skin depth and rigorously derive a spectral envelope from the model via asymptotic 
analysis. 

Especially for p-polarized oblique incidence, the formation of highly dense and narrow electron nanobunches in 
front of the surface is often observed. If these bunches carry a considerable amount of charge, they emit intense 
high frequency radiation that is not described within the ROM model. In this case, we can derive the spectrum by 
calculating the coherent synchrotron emission (CSE) from these bunches - as presented in subsection III D 



A. Starting point of analysis 



The foundation of the theory of surface HHG is explained here. 

In Sub. |III A 1[ we demonstrate the Green function solution of the inhomogeneous wave equation 
provides a general starting point of HHG theory. Then (Sub 



This solution 

HI A 2| , we have a closer look at the source term to 
derive general selection rules concerning the parity and polarization of the generated harmonics. In Sub. |III A3[ we 
briefly deal with surface HHG in the sub-relativistic regime. Sub. |III A 4| explains the meaning of the ultra-relativistic 
S'-similarity group for HHG. 



Solution of the inhomogeneous wave equation 



Let us begin with the classical wave equation for the electromagnetic potential. Since the basic structure of the 
physical mechanism can best be understood in a one dimensional slab geometry, we let A depend on only one spatial 
coordinate x and the time t. By this we do not generally exclude oblique incidence, because it can be treated in a 
Lorentz boosted frame wherein the laser is normally incident (see App. [B]). In this geometry, the wave equation in 
Coulomb gauge (V • A = 0) can be written as: 



^d^A{t,x)-dlA{t,x) = —u{t,x). 
c' c 



(1) 



Equation ([T|) can generally be solved with the help of a Green function. We formally write down the solution as 

A{t,x)=4Tr if jj_{t',x')G{t,x,t',x')dt'dx'. (2) 



By the choice of G, the asymptotic behaviour of A can be controlled. We choose the Green function G in order to 
solve Eq. ([T]) and additionally fulfll the boundary condition |A(f, a;)| for x +oo, i.e. there is no light coming 
from the right and all radiation coming from the left is fully reflected. We obtain: 



G{t,x,t',x') = 



t-t' 



\x~x'\ 



-e[t-t' - 



(3) 



where denotes the Heaviside step function. 

It proves convenient to continue working with the transverse electric fleld E_l — —c'^dtA instead of the vector 
potential here. Thus Eq. ^ becomes: 



^±{t,x) 



2tt 



j_lU ,x -J_L tH ,x 



dx'. 



(4) 
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Let us define Xy = sup {x : j{t, x') = 0, \/t, Vx' < x}: the leftmost point which is reached by any charge during the 
laser-plasma interaction process. It can be seen that for all a; < to the left of the plasma, the first term in Eq. Q 
represents the incoming radiation, while the second term represents the refiected one. To the right of the plasma both 
terms cancel, as our choice of the Green function requested. 

Due to the assumption of one-dimensionality, the radiation does not change while propagating in vacuum, and the 
incoming and outgoing fields and are each function of only one variable t ± x/c. We may therefore drop the 
argument x and identify 

E,(t) = — / j^[t+-,x'] dx' (5) 
E,(i) = / j±[t--,x'] dx', (6) 



so that 'E±{t, x) —'Ei{t~ x/c) +Er{t + x/c) for x < Xy. 

Eqs. ([5| and ^ provide a powerful and general starting point for the theory of harmonics generation. Whereas 
Eq. (|6| tells us how to obtain the refiected field E^ from a given current distribution, Eq. ([s]) provides a condition on 
the current for a given incident field E^ . Note that this equation stems from our choice of the Green function and 
physically represents the condition of total refiection. If instead we had chosen the Green function in a way that all 
fields vanish for t — — oo, then there would be no E^, but the generated field to the left of the plasma would be the 
same as E^ in Eq. Q. Equation ^ is of course not sufficient to explicitly calculate j_L, but it can be harnessed to 
obtain j and consequently Ej, in conjunction with some additional assumption. This is a possible approach to derive 



boundary conditions for the ROM and TROM models (subsections III B and III C I 



2. Selection Rules 



Before we move on to present these models, let us collect some general facts about the source term jj^, stemming 
from the plasma response to the laser. Therefore, we consider the ffuid equations for a cold relativistic plasma. These 
equations do not account for kinetic effects like trajectory crossing, but they are adequate to derive some general 
properties of the physical process. For a relatively short laser pulse, we can neglect the ion response, so the current 
is given by: 

j = {nv - novo) , (7) 

where e is the elementary charge, n is the electron density and v is the electron ffuid velocity. Note that we consider 
all magnitudes in the inertial frame in which the laser is normally incident. In this frame, the electrons and ions 
possess some initial velocity Vg parallel to the surface, and the initial density uq is not necessarily identical to the 
initial density in the laboratory frame. The velocity v is related to the relativistic momentum p like v = p/jnie, 
where 7 = \/l + (p/toec)^- Due to the conservation of the canonical momentum|23], the transverse component can 
directly be connected to A in the presumed ID geometry: 

P_L = Po + eA. (8) 

The set of equations is completed by the equation of motion for the longitudinal momentum component p^, the 
continuity equation and the Poisson equation for the electrostatic potential due to charge separation: 

dtn = {nv^) , (10) 
d^<P = 47re(n-7io), (11) 

wherein d/dt — dt + v^dx denotes the absolute time derivative. 

Having a closer look at these equations, it is possible to derive some "selection rules" with respect to parity (even 
or odd harmonic numbers) and polarization (see also Ref. [15 ). 

1. Assuming normal incidence of linearly polarized light, we take po = 0, A = Aey. It is obvious then, that jz — 0, 
and the polarization of the incident light is conserved. For the y-component of the source term, we obtain 



jy oc nA/^. About the longitudinal momentum p^, which enters both n [through in Eq. (10l[ and 7, we 
know that Px = e [dx^ + cA^ /j^j. Thus, the longitudinal momentum is driven by the square of the laser ffeld 
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incident light 


odd harmonics 


even harmonics 


normal (linear) 


same as incident 




oblique (s) 


s 


P 


oblique (p) 


P 


P 



Table I: Selection rules for polarization (s, p) and parity (even, odd) of harmonics at plasma surfaces depending on the 
polarization and the angle of the incident laser. 



and therefore has the same periodicity as A^. Consequently, it possesses only even harmonic s of the laser 
frequency. The same holds true for 7i and 7, as can be seen from Eq. ( 10 1 and 7 = -y/l + eA^ + p'^. Finally, jy is 



a product of A ^ cosojo^ (in zeroth order) with quantities that possess only even harmonics of the fundamental 
laser frequency. We conclude, that jy and therefore A purely consist of odd harmonics of the fundamental. 

2. For s-polarized oblique incidence, we may assume po = poe^ and initially A — Aey. In this case, jy cx nA/j 
as in the normal incidence case, but additionally there is a source term in z-direction: jz cx npo/j — noPo/7o- 
Again, p^., n and 7 contain only even harmonics of the laser frequency. Consequently, jy and Ay contain only 
odd harmonics and jz and Az contain only even harmonics of the fundamental. 

3. For p-polarized oblique incidence, we can take po = po^y and A = Aey. We immediately see, that there is no 
source term in z-direction (jz = 0), and jy oc pq {n/"f — no/70) + enA/j obviously contains both even and odd 
harmonics. Another interesting fact is that p^ is now also driven by a term that is linear in A. This implies, 
that harmonics can be observed here at lower intensities compared to s-polarized and normal incidence. 

In table [T] the rules just derived are summarized for reference. 



3. Sub-relativistic plasma non-linearity 

This work deals with generation of harmonics due to relativistic mechanisms. These have to be distinguished from 
harmonics generated by sub-relativistic plasma non-linearity. Here, we explain the sub-relativistic mechanism in brief. 

It is found that for p-polarized oblique laser incidence, the threshold for harmonics generation is much lower than 
for s-polarized or normal incidence. This is due to plasma non-linearities, which are not of relativistic origin and only 
occur for p-polarized incidence. Under this condition, two effects may lead to the excitation of plasma oscillations 
inside the inhomogeneous plasma-gradient: 

1. Resonant absorption of the laser field, see e.g. the book by Kruer [23]. 

2. Electron bunches that are separated from the main plasma and then re-enter, see the famous work by Brunei 

Due to the strong inhomogeneity of the plasma, these oscillations couple back to electromagnetic modes via sum 
frequency generation, leading to the emission of high harmonics. When the excitation happens by means of Brunei 
electrons, the mechanism is commonly referred to as "coherent wake emission" (CWE) [19, 24_. CWE is the prevalent 
sub-relativistic generation process for femtosecond-scale laser pulses. 

According to their generation mechanism, the sub-relativistic harmonics have a strict frequency limit, given by 
the plasma frequency Wp corresponding to the maximum density [TJ IHl [12] • The subsequently discussed relativistic 
harmonics are not subject to this limitation and can therefore easily be distinguished from the ones generated by the 
non-relativistic mechanism. The transition between both regimes for moderately relativistic laser pulses was discussed 
by Tarasevitch et al. in Ref. |25| . 



4-. S -similarity 

The 5-similarity ^264 is a similarity group that characterizes the interaction of higly relativistic (oo ^ 1) lasers with 
plasmas. The major dimensionless parameter here the similarity number S = Ne/aoNc- If S is kept constant, but ao 
and TVe are changed, the behaviour is similar, so that many important physical quantities can be obtained by simple 
scaling laws. The theory was first applied to laser - underdense plasma interaction, where it delivered useful scaling 
laws for laser- Wakefield acceleration (LWFA). 
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Concerning the interaction of lasers with overdense plasmas, similarity theory has to be applied with care, since 
in the skin layer with its extreme field gradients, ultra-relativistic and weakly or non-relativistic electron motion 
happens in directly adjacent regions. Even if the laser is highly relativistic in the sense oq ^ 1, most of the electrons 
deeper inside the skin layer move with only moderately relativistic velocities as long as 5 > 1. For S < 1, relativistic 
transparency sets in. Nevertheless, some signatures of S'-similarity can be observed here. 

According to similarity theory, the momenta of the highly relativistic electrons can be written as 

p^iao,N^;t)^aoMS;t), (12) 

where is a characteristic function that describes the motion of the plasma electrons. Without knowing the details, 
we can assume |pi| ^ 1 close to the surface, i.e. a considerable fraction of the laser field accelerates the electrons at 
the surface. This is certainly the case in parameter ranges where harmonics are efficiently generated. 

Let us now consider the x-component of the velocity. Dropping the particle index i for simplicity, we write = 
Px(ao ^ +pI + pD^^^^ and consequently 

- v/T^ V «o-^ + pi 

_ I aoPj; (when p^ = 0) ^^^^ 
1 ^Jl + p1/i>\ (otherwise) 

It is evident that grows very large - up to the order of oq - around the times, when the transverse momentum 
components vanish. Further, we see that must possess a distinct maximum here - a kind of spike - since otherwise 
^x remains in the order of one and oq 3> 1. Later we show, that the generation of high order harmonics happens 
basically due to these "7-spikes". 

Note, that this also implies that the generation of relativistic harmonics is generally much more efficient for p- 
polarized and normally incident light than it is for s-polarized light. In the case of s-polarized oblique incidence, the 
momentum space is three-dimensional and there are two generally non-vanishing transverse momentum components. 
Consequently, the trajectory of p in momentum space does not necessarily cross the Pa,-axis and 73. probably remains 
in the order of 1 during the entire interaction process. For p-polarized or normally incident light, one of the transverse 
momentum components vanishes for symmetry reasons, = 0. Therefore, the momentum space is two-dimensional 
and p has to cross the p^j-axis every time py changes sign. 



B. The relativistically oscillating mirror (ROM) model 



Due to its descriptive nature, the term "relativistically oscillating mirror" (ROM) is in common use. However, its 
usage varies among authors and there has been no accurate and generally accepted definition so far. In the frame of 
this work, we define the ROM model as the model based on the boundary condition 



E,, t 



XAKp{t) 



Er t 



= 0, 



(14) 



wherein xarp denotes the coordinate of the "apparent reflection point" (ARP). In Ref. [15 , where the term "oscillating 
mirror" was first used in the context of relativistic laser-plasma interaction, it was applied to a model based on an 
oscillating step-like plasma boundary. We are going to see soon that the above boundary condition is closely related to 
that model. Further, the ARP is intuitively understood as a sort of mirror, which oscillates at relativistic velocities. 

This subsection consists of three parts. At first (Sub. IIIB 1 ), we investigate the foundation of Eq. ( 14l, trying to 
clarify, under which conditions it is applicable. Then (Sub. IIIB 2 1, we demonstrate a simple way to check the validity 
of the model within a simulation. Finally (Sub. IIIB 3), we derive some very general properties of the spectrum that 



follows from Eq. ( 14 1 



Foundatton of the ARP boundary condition 



Here, the applicability of the boundary condition Eq. 
to arrive at the condition. The first approach was introduced by Gordienko 
Taylor expansion of the current distribution. The alternative approach connects the ARP boundary condition to the 



(14 1 is analyzed. To do this, we consider two possible ways 
et al. 126 in 2004. It is based on the 
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assumption of a moving step- like electron density profile. The correlation between the shape of the electron density 
profile and the resulting radiation can be confirmed within PIC simulations. 

We begin with the approach from Ref . [51 US] • It is based on the Taylor expansion of the current distribution with 
respect to time: j_L(io + h, x) « j_L(to, + hdtjj_{to, x) + 0{h?). This is inserted into Eq. Q. Expanding around 
to — t — t' and using h — {x — x')/c, we find that the zeroth order term vanishes immediately. Keeping the first order 
terms and neglecting the second and higher orders we get 

E^^— ^^dti^it - 1', x') dx'. (15) 

C Jx C 

Now, the time derivative dti± is estimated by the current divided by the "skin layer evolution time" r = 
min (j_L/9tjj^). In this way, for x ~ Xgurf at the plasma surface, Ej^ can be approximated as: 

47r(5 

E_L(a:: = Xsurf) ^ J_L, (16) 

CT 

where b refers to the skin length and Jj^to the instantaneous net current. The characteristic time r can be estimated 
by the inverse laser pulse frequency, r ^ <^(7^- In the linear approximation, the skin length is given by the plasma 
frequency, b = c/u)p. Thus, for highly overdense plasmas, we expect the field being small at the plasma surface and 



may consequently apply Eq. ( 14 1 . 

In order to further clarify the conditions which lead to the fulfilment of Eq. ([l4 1 , one can also consider an alternative 
derivation of the boundary condition (14). This derivation is based on the model of an oscillating, step-like boundary. 
We start with an arbitrary polarization component of the wave equation ([T|). The equation is then adapted to the 
step-like density profile and normal incidence. Further, we make use of the canonical momentum conservation ([8| and 
switch to relativistically normalized units {ct — t, ...) for convenience: 



w2 



dt -^dt\A = e{x- x,urf(i)) -^A, (17) 
where ujp is the electron plasma frequency and 7 is the electron 7-factor. We make the complex ansatz: 

A{t x)^\ ^'''^ ~ ^/^-^ ^ ^''^^ ^ ^^^"^ a;surf(t)) (-j^g^ 

yAs{t + iKxlc) {X > Xsnrf{t)) 

wherein k = ^Ji-^p/ (71^0 ) ~ ^ ^ '^^^l number, as the plasma is overdense. To take account for the relativistic non- 
linearities, we allow for general functions instead of strictly assuming Ar, A^ (x exp [iuiQt). Note, that the vacuum part 



{x < Xsurf) of Eq. (18 1 is an exact solution of Eq. (17). The skin layer part is an exact solution for the fundamental 
mode. As oc exp (iujot). Taking into account that in many cases (w~^/^ - spectrum) the biggest share of energy is still 
contained in the laser fundamental mode, we consider this approximation being reasonable enough within our simple 
model. 

Now the function A as well as its first partial spatial derivative must be continuous at the point a;surf(0 at every 
time t. Defining ai{t) = eAi(t — XsuYfit)) /mc, ar{t) = eAr(t + Xsmf(i)) /mc, as{i) = eAs{t -\- iKXsuY{it))mc, we get: 

ai + ar — fls (19) 
111 

"flr — 1 71— Is- (20) 



Solving for yields 



isurf/c - 1 isurf/c + 1 Xsurf/c - i/k 



as = E,it x,Ut)/c), (21) 

1 + 1/ K 



where di = (isurf/c — I) Ei was used. Inserting this back into Eq. (20), we obtain 



Er{t + x,^rf{t)/c) + ^-^E,{t - XsuAt)/c) = 0. (22) 

K + t 
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Figure 1: Radiation in time (a) and spectral (b) domain in tlie ROM regime. Simulation parameters are: normal incidence, 
plasma density A^^ = 250 Nc, laser field amplitude is ao = 60. In both frames, the reflected field is represented by a blue line. 
In (a), the green line represents the field of the incident laser and the black dashed lines mark the maximum field of it. In (b), 
the dotted black line represents an 8/3 power law. 



Since \{k — i) / {k + i)\ = 1 for k e M, it is now seen that Eq. (22 1 agrees with Eq. (14 1 except for a phase term, 
phase can be included in the function XARp(t), setting 



1 

UJQ 



arccos 1 — 



This 



(23) 



We have shown here that the ARP boundary condition ( 14 1 is valid under three main assumptions: a step-like 



electron density profile, normal incidence and the interaction is dominated by the laser fundamental. In this case, the 
ansatz (18) is reasonable. If the density inside the plasma is not exactly constant, but there is a sharp rising edge 
behind which comparatively weak fluctuations follow, the ansatz ( 18 1 might still be useful, as the precise behaviour 
of the fleld deep inside the skin layer has no strong influence on the reflection. The new calculation equips us with a 
rough idea of when Eq. ( 14) can be expected to be useful. 



2. Verifying the ARP boundary condition 



Equation ( 14 1 has a simple interpretation that allows us to verify within simulation data, whether it is fulflUed or 
not. For Eq. (14) to have any useful physical meaning, it is required that |iARp(i)| < c at all time. Otherwise, the 
fields would be un- or overdetermined, leading to contradictions or useless tautologies. 

It is obvious then, that the reflected fleld is nothing but a phase modulation of the negative of the incident 
one {—Ei). In a PIC simulation, we can easily check this by looking at the fields in the time domain. If and only if 
Eq. ( [l4|) is fulfilled, then both functions possess the same sequence of extrema and monotonic intervals. 

Fig. llTa) shows an example, where Eq. (14) is fulfilled to a good approximation. This example was obtained by 
normal incidence on a sharply defined plasma, therefore confirming the considerations presented above. However, 



compare this to Fig. |6ja) to see that Eq. (14) cannot always be applied. In Sec. HID we are going to discuss this 
interesting case in more detail. Now we continue with the case, where the ROM model works. 



3. Analytical derivation of the spectrum 



Given the validity Eq. ( 14 1 , it is possible to calculate the general form of the spectral envelope with only a few 
straightforward assumptions. 



We start by writing down the Fourier transform of from Eq. (14), taking into account the retarded time: 



Erioo) = - E,[t- ^) e-(*+-/^) 1 + 



dt. 



(24) 
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The incoming laser pulse is described by an envelope approximation 
/(t) [exp {i{u;ot + 0o)) ~ exp {—i{ujQt + (po)] /2. where g(t) is a slowly varying function. We arrive at 



E^it) 



Eriuj) 



E± = - / .9 - 



XARP 



C 



X exp 



t \ LU [t -\ ] ±LUq [t 



(25) 

Now note that for high u, the exponential term leads to a rapid oscillation of the integrand during most of the 
time. Because of this oscillation, most contributions cancel, except for those where the phase of the integrand stands 
still. This means that the integral can be handled by the method of stationary phase. The somewhat more technical 
details of this calculation are shifted to appendix]^ but before presenting the final result, we would like to remark 
two interesting points: 

First, the stationary phase points correspond to the instants when the ARP moves towards the observer with 
maximum velocity. These moments are crucial for the generation of high order harmonics. The corresponding ARP 



1 



^ARP 



-1/2 



possesses a sharp spike at these instants, which is the reason why we also 



gamma factor 7arp = 
call them j-spikes 

Second, the spectrum depends on the exact behaviour of the ARP in the neighbourhood of these points. In Ref. jB] 
it was presumed, that the derivative of the ARP acceleration is different from zero at the 7-spike. Other cases are 
imaginable however, and it is intriguing to see, what difference they make. We here consider the most general case, 
in which xarp has a zero of order 2n — 1 at the 7-spike, meaning that d'^XARp /dt*' — for all 2 < A: < 2n. We will 
subsequently refer to n as the order of the ^-spike. 

After the calculations in appendix \K\ the spectrum can be written as: 



■1 + 4 
n + 1 



'^£{-1,1} 



'gAi„ 



a;7 







,2 {auj) 



l/(2n-fl) 



(26) 



wherein 7 refers to the peak value of 7arp(0 and a is a constant related to the behaviour of the ARP trajectory 
close to the 7-spike. As the ARP motion happens on the timescale of the laser period, we may estimate a 



2n 



UJQ 



and consequently write a = aojQ , where a is a numeric constant on the order of one. (po is connected to the phase, 
at which the 7-spikes occur. As we are going to see later, it has no signifcant influence on the spectrum, unless 
(/)o ~ (n + l/2)7r. gAi„ is a generalized Airy-function as defined in appendix [C| These functions are not commonly 
available in general purpose numerical function libraries. With a small trick, they are however not hard to compute. 
The details of the numerical calculation of the integral are explained in App. [C] We now regard the common case 
n = 1 and the special case n > 1, corresponding to instantaneously vanishing acceleration at the 7-spike, separately. 
The most typical case n = 1 has been investigated by Baeva, Gordienko and Pukhov (BGP) in Ref. [6J. In this case. 



Eq. (26) can be written with a conventional Airy function: 



^(w) 



Ai 



2 (aw) 



1/3 



(27) 



The spectrum is presented in Fig. [2] We see that they consist of a power law which "rolls over" into an exponential 
decay at a frequency, which depends on 7. By comparing the argument of the Airy function to one, we find this 
characteristic roll-over frequency to be 



UJr — vSaLdQ'J 



(28) 



Note the favourable scaling of this characteristic frequency cj^ with 7^^, in contrast to the Doppler shift from the 
refiection at a constantly moving mirror, which produces a frequency upshift by a factor of only 47^. 

Below this frequency, the airy functions are almost constant and the spectrum can be approximated by the famous 
power law 



BGP 



(w < UJr) 



1 



(29) 
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Figure 2: The BGP-spectmm Eq. (27 1 in a double logaritlimic representation. The dashed black line denotes an 8/3-power law 
and the dotted line represents the roll-over frequency ujr = \/8(5a;o7^- For the picture, we took a = 1, 7 = 5 and (t>o — 0, but 
the resulting spectra are largely independent of these variables, provided the cj-axis is scaled appropriately. 




Figure 3: Experimental measurement of the relativistic harmonics spectrum by Dromey et al., taken from Ref. 20]. The dots 
indicate measured harmonics normalized to the 238th harmonic, the lines correspond to power law fits / oc with q = 2.5 
for the best fit (red line), and q — 2.2 (blue) respectively q = 2.7 (black), close to the theoretical value q — 8/3 of the BGP 
spectrum Eq. (271. 



This scaling can nicely be seen in Fig. [2] Technically, there is an exception to this scaling: the special case 
(pQ Ki [n + 1/2)tt. In that case, the leading orders of the Airy functions cancel each other out, yielding a slightly 
steeper power law decay: /bgp ^ w^^^^'^, see also Fig. 5[b). This special case seems to be of not much physical 



importance however. 

For frequencies much larger than ujr, the Airy function dominates and the decay becomes exponential: 

/BGp(w>Wr)- f — ) '^CXpf-^— ) . (30) 

An example for this sort of spectrum from a PIC simulation can be seen in Fig. [ijb). An experimental spectrum, 

measured over a wide range of frequencies, is shown in Fig. [3] It is found to be a power law, and the exponent 

- - - ' — ■ - - - - - 3 



q G [2.2; 2.7] agrees with the BGP spectrum (27 1. Further, the scaling of the cut-off frequency was found to be ~ 7 



also in agreement with the theoretical predictions of the discussed model. 
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Figure 5: Logarithmic plot of the spectra Eq. (26 1 following from the ARP boundary condition, for orders n = 1, 2, 10 of the 
7-spike. (a) shows the case (jio ~ 0, (b) shows the special case (f>o ~ vr/2. Unless (po ~ {n+ l/2)7r, the observed spectra resemble 
(a) rather than (b). All spectra have been normalized to 7o = /(i^o), and a = 1, 7 = 8 was used throughout. 



Let us now move on to the cases 71 > 1. Fig. [5] shows the spectra Eq. (26 1 for different orders n of the 7-spike. To the 



best of our knowledge, the cases n > 1 have not been investigated before. 

In Fig. [Sj^a), the typical case (f>o — is displayed. Fig. [sf^b) shows the case 0o — l^i where the first order term of 
the generalized Airy-functions just cancel out. It makes sense to also consider this special case here, as there is reason 
to believe that the occurence of higher-order 7-spikes may be related to the occurence of the spikes at a certain phase. 
Comparing the spectra for the higher order 7-spikes n > 1 to the BGP case we note the following differences: 

1. The power law part of the spectra decays slightly slower with increasing n. In addition, particularly in the case 
(/•o = ■''■/2, the non-power law part gains influence for increasing n at low frequencies already and leads to even 
slower decaying spectra. This is possibly favourable for the efficient production of attosecond pulses. 

2. Because of the oscillatory behaviour of the generalized Airy functions gAi„(x) at positive x and for n > 1, the 
spectra become strongly modulated at frequencies u) > uir-, compare Fig. [5] Also, this effect is more pronounced 
for (/)o — "K 12. This might explain some of spectral modulations observed in numerical and real experiments 
before, see e.g. Ref. [57]. 

3. The roll-off frequency, which scales as 7"^ in the BGP case, approaches a 7^ scaling in the limit n — ^ cx3, 
reminiscent of the Doppler effect from a mirror moving with constant velocity. This seems reasonable, since 
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for higher order 7-spikes, the acceleration is very small in the neighbourhood of the stationary phase point. 
Therefore, its influence on the spectrum decreases with n. 

Using the estimates from Sec . |III A4] based on the S'-similarity theory, we can further assume that for fixed S'-parameter 
7 oc Uq. In the BGP-case n = 1 this means oc Oq, which is in agreement with experimental observations reported 
in Ref. [20,. 

To sum up this subsection, we have reviewed the popular "relativistically oscillating mirror" (ROM) model for 



the relativistic generation of harmonics at overdense plasma surfaces, based on Eq. (14 1. We have found, that it is 
applicable for normal incidence and step-like electron plasma boundaries. Because of its simplicity - it reduces the 
whole complex interaction physics to one simple function XARp{t) - the model helps us to gain insight into the basic 
mechanism that leads to the generation of high harmonics. Further, we have analytically calculated the possible 
spectra in the relativistic limit with the help of asymptotic analysis. Here we noticed, that even within the model, 
spectra that deviate considerably from the well known BGP 8/3-power law are in principle possible. 



C. Totally reflecting oscillating mirror (TROM) and 



-2/3 



spectrum 



As we see e.g. from Fig. [6] the ROM model based on Eq. (14 1 as it was used in Ref. [6l, is not always valid in the 
highly relativistic regime - not even as an approximation. It is thus worth looking for alternatives. 

This subsection is about another model that one might intuitively associate with the name "relativistically oscillating 
mirror", we call it the "totally reflecting oscillating mirror" (TROM). The model is rigorously based on the assumption 
of total reflection from a perfectly localized current layer. 



First (Sub. Ill C 1 1, we derive the corresponding boundary condition. Then (Sub. Ill C 2 ), we demonstrate the spec- 
tral properties that follow from this boundary condition, utilizing asymptotic analysis once again. Finally (Sub. Ill C 3 1 
we give some remarks about the possible physical realization of the model. 



1. Foundation of the TROM boundary condition 

The TROM model is particularly interesting because of its mathematical lucidity. It can be rigorously derived from 
only two straightforward assumptions. These assumptions are: 

1. There is total reflection, no light passes through the mirror. Therefore, we can relate the plasma current to the 
incident radiation via Eq. ([s]). 

2. The skin layer of the reflecting plasma is infinitely thin. Therefore, the current can completely be described by 

= j{t)6{x - XTnoM{t)). 

Inserting the current profile into Eqs. ([5| and (|6|, we obtain: 

m) = -2vr (31) 

EAt) = 2. J(^ + 7hom(0/c) (32) 
c - a;TROM(i + a;TROM(i)/c) 

Now, the assumption of total reflection is exploited by using Eq. Q. We eliminate j and arrive at the boundary 
condition: 

eJi^ ^^'^"^^^^ \ + ^- E,(t- 

\ c / l + a;TROM/c V c y 

Compare this to Eq. ( [l4| ), which represents the ROM model. The difference lies in the pre-factor of Ei, which 
amplifies the refiected field at times, when the mirror moves towards the observer. Since these are the regions which 
are responsible for high frequency radiation, we expect a flatter spectrum here compared to the ROM model. 



Further note, that Eq. (33 1 is always the correct boundary condition for a totally reflecting mirror in the limit of 



constant velocity. In this case, Eq. (33 1 could simply be derived by a Lorentz transformation to the system, where 



the mirror is at rest. For a strongly accelerated mirror however, we need the additional assumption of a perfectly 



localized skin layer to obtain Eq. ( 33 1 
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2. Analytical derivation of the TROM spectrum 



We now derive the spectrum corresponding to Eq. (33). The beginning of the calculation is analogue to the 
calculation in subsection |III B 3| and we arrive at 



^4 



a:TROM \ 
t J exp 



t + 



c J 



a^TROM 



(34) 



Compare this to Eq. ( 25 1 . The difference lies in the last factor: Whereas in Eq. ( 25 1 it is 1 + xarp /c, we have 



1 — a;TROM/c here. This difference is crucial, since at the stationary phase points, where iTROM/ARP ~ ~c, the term 
in Eq. (251 becomes very small, whereas the term in Eq. (34 1 even has a maximum. 



Again, we can analytically calculate the corresponding spectrum, as shown in App. In general, for a 7-spike of 
the order n we obtain: 



^trom(w) oc 



J2 '^S^'^r 



W7 



a4:(jj 







<Te{-i,i} 



l/(2n+l) 



(35) 



This is the same as the ROM spectrum Eq. (26 1, except for the different exponent in the power law. The TROM 
spectrum is much flatter. For high order 7-spikes, the power law part even tends to a;", so that the spectrum is merely 
determined by the generalized Airy functions. 

In the more likely case n = 1, Eq. (351 can to a good approximation be simplified: 



-^trom('^) 



,2/3 



Ai 



2/3\ 



(36) 



Compared to the uj decay predicted for the ROM model [Eq. (27l], we obtain a slowly decaying uj ^/'^ power 
law here. 



3. Physical Feasibility of the TROM model 

As we have seen, the TROM model yields a distinctly flatter spectrum than the ROM one. Therefore, if there were 
a physical system that behaves according to the TROM model, it could be much more efficient in the production of 
attosecond pulses. Let us try to answer (a) why this is difficult and (b) how it might still be possible. 



The difficulty can readily be seen from Eq. (31). We notice that the current j does not necessarily vanish at 
the instant when the surface moves at maximum velocity. This is in contrast to the normal behaviour of an ultra- 
relativistic plasma. The transverse current is the product of the transverse fluid velocity component Vy and the 
charge density p. Since the transverse velocity component becomes very small at the instant of maximum longitudinal 
velocity, a flnite j implies a huge plasma density. But very dense plasmas are hard to drive to relativistic motion. 

For single pulse schemes, the realization is probably impossible. The behaviour of ultra-relativistic plasmas is 
governed by the S'-parameter S = Ng/aoNc- If the S'-parameter is too low, it leads to an extended skin layer in 
contradiction to the assumption of a perfectly localized current layer. If the S'-parameter is too high, the plasma is 
not driven to relativistic motion at all. 



In Ref. |22| , Tarasevitch et al. propose the realization of the boundary condition (33) via a two pulse scheme. In 
the scheme, the flrst, relativistically strong pulse drives the plasma surface to oscillation. The second pulse is much 
weaker and has a polarization orthogonal to the flrst one. It is used as a probe and the spectrum in the direction 
of its polarization is recorded. Indeed, for a certain set of parameters it was possible to observe the generation of 



harmonics according to Eq. (33). Thereby, they heuristically also flnd a 2/3-power law spectrum, conflrming the 
above calculations. 

In this case, the probe pulse "harvests" the harmonics generated by the much stronger driver pulse. Thus, the 
scheme is not appropriate to increase the overall efficiency of frequency conversion or attosecond pulse production. In 
the following section, we are going to look at a physical mechanism, where the overall efficiency is indeed increased 
considerably in comparison to the ROM case. 
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Figure 6: Radiation in time (a) and spectral (b) domain for simulation in the "nanobunching" regime. Compare this to Fig. [T] 
Simulation parameters are: plasma density ramp oc exp(a;/(0.33 A)) up to a maximum density of — 95 Nc (lab frame), 
oblique incidence at 63° angle (p-polarized) , ao — 60. Again, the reflected field is represented by a blue line, the green line 
represents the field of the incident laser and the black dashed lines mark the maximum field of it. In (b), the dotted black line 
represents an 8/3 power law and the red dashed line corresponds to the analytical nanobunch CSE spectrum given by Eqs. ( |41[ ) 
and (421, with ujrs = 800 cuo and cOrf = 225 wo . 



D. Coherent synchrotron emission (CSE) from electron nanobunches 



Cases where the ARP boundary condition (14 1 does not apply are studied here. We find, that the radiation can 



be described as coherent synchrotron emission (CSE) from extremely compressed electron "nanobunches" that form 
in front of the surface. 



At first (Sub. HID 1 1, the generation process is investigated by close examination of PIC data. Then, the spectrum 

Finally, we analyze the sensitivity of the process to changes in the laser-plasma 



is calculated analytically (Sub. Ill D 2 ) 
parameters (Sub. HID 3). 



1. Electron nanobunching process 

Let us now have a fresh look at Fig. [6j It is evident, that the maximum of the reflected field reaches out about 
an order of magnitude higher than the amplitude of the incident laser. The reflected radiation can clearly not be 



obtained from the incident one by phase modulation and the ARP boundary condition Eq. ( 14 1 fails. 

Consequently, the spectrum deviates from the 8/3-power law, compare Fig.|6jb). Indeed, the efficiency of harmonics 
generation is much higher than estimated by the calculations in Ref. [6]: about two orders of magnitude at the 
hundredth harmonic. Also, we can securely exclude coherent wake emission (CWE) as the responsible mechanism, 
since this would request a cut-off around uj = lOwo- The radiation has to be attributed to a new sort of mechanism. 

To get a picture of the physics behind, let us have a look at the motion of the plasma electrons that generate the 
radiation. Figure [7] shows the evolution of the electron density corresponding to both sample cases from Figs. [T] and [6] 
In addition to the density, contour lines of the spectrally flltered reflected radiation are plotted. These lines illustrate 
where the main part of the high frequency radiation emerges. 

We observe that in both cases the main part of the harmonics is generated at the point, when the electrons 
move towards the observer. This shows again that in both cases the radiation does not stem from CWE. For CWE 
harmonics, the radiation is generated inside the plasma, at the instant when the Brunei electrons re-enter the plasma 

m- 

Apart from that mutuality, the two presented cases are very different. Figure ma) corresponds to the ROM case. 
It can be seen that the density proflle remains roughly step-like during the whole interaction process and the plasma 
skin layer radiates as a whole. This explains why the ROM model works well here, as we have seen before in Fig. [T] 

Figure |7|^b) looks clearly different. The density distribution at the moment of harmonics generation is far from 
being step-like, but possesses a highly dense (up to ^ 10000 A^c density) and very narrow (5-like peak, with a width of 
only a few nanometres. This electron "nanobunch" emits synchrotron radiation coherently. 
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Figure 7: The electron density and contour lines (cyan) of the emitted harmonics radiation for lj/ixjo > 4.5, in (a) the ROM 
and (b) the nanobunching regime. The small windows inside the main figures show the detailed density profile at the instant of 
harmonic generation. All magnitudes are taken in the simulation frame. The simulation parameters in panel (a) are the same 
as those in Fig. [l] and (b) corresponds to Fig. [g] 

The high frequency radiation is emitted by a highly compressed electron bunch moving away from the plasma. 
However, the electrons first become compressed by the relativistic ponderomotive force of the laser that is directed 
into the plasma, compare the blue lines in Fig. [8] During that phase, the longitudinal electric field component grows 
until the electrostatic force turns around the bunch, compare the green lines in Fig. |8] Normally, the bunch will loose 
its compression in that instant, but in some cases, as in the one considered here, the fields and the bunch current 
match in a way that the bunch maintains or even increases its compression. The final stage is depicted by the red 
lines in Fig. |8] 

We emphasize, that such extreme nanobunching does not occur in every case of p-polarized oblique incidence of 
a highly relativistic laser on an overdense plasma surface. On the contrary, it turns out that the process is highly 
sensitive to changes in the plasma density profile, laser pulse amplitude, pulse duration, angle of incidence and even 
the carrier envelope phase of the laser. For a longer pulse, we may even observe the case, that nanobunching is present 
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Figure 8: Formation of the nanobunch in the simulation corresponding to Fig. [6]and[7|b). We depict the electron density Ne 
in units of the critical density Nc, the transverse magnetic field component and the longitudinal electric field component 
Ex in relativistically normalized units. 

in some optical cycles but not in others. The parameters in the example were selected in a way to demonstrate the 
new effect unambiguously, i.e. the nanobunch is well formed and emits a spectrum that clearly differs from the BGP 
one. The dependence of the effect on some parameters is discussed in subsection |f 1 1 D 3l 

Because of the one dimensional slab geometry, the spectrum is not the same as the well known synchrotron spectrum 
[5S] of a point particle. We now calculate the spectrum analytically. 

2. Analytical derivation of the nanobunch ID CSE spectrum 
The calculation of the spectrum is based on two assumptions: 

1. As in the TROM model, the radiation is generated by a narrow bunch of electrons. Optimal coherency for high 
frequencies will certainly be achieved, if the current layer is infinitely narrow: j{t,x) = j{t)5{x — Xei{t)). To 
include more realistic cases, we allow in our calculations for a narrow, but finite electron distribution: 

j{t,x)=j{t)f{x-x,i{t)) (37) 

with variable current jit) and position Xei{t), but fixed shape f{x). 

2. In contrast to the TROM model, we give up on calculating j(t) directly from the incident radiation by the 
assumption of total reflection. This means, that although we know that the bunch itself is not capable of totally 



17 



reflecting the incoming radiation and consequently there are some additional currents inside the plasma, we do 
not care for them as their contribution to the high frequency spectrum are small compared to the contribution 
by the highly compressed bunch. 

However, to get some kind of result, an assumption about the functions j{t) is required. Since we are dealing 
with the ultrarelativistic regime ao S> 1, it is reasonable to assume that changes in the velocity components 
are governed by changes in the direction of motion rather than by changes in the absolute velocity, which is 
constantly very close to the speed of light c. We are going to see, that this assumption is enough to obtain the 
spectrum. 

Following Eq. ([6|, the radiation field is expressed as i?csE a^) — 27rc^^/ j (t + (x — x')/c, x') dx'. We take the 
Fourier transform, thereby considering the retarded time, and arrive at the integral 



2tt - 
c 



+00 



j(t) exp 



'ioj t + 



dt, 



(38) 



wherein f{uj) denotes the Fourier transform of the shape function. 

In analogy to the standard synchrotron radiation by a point particle, the integral can be solved with the method of 
stationary phase. Therefore, we Taylor expand the current j{t) and the electron bunch coordinate Xei{t) around the 
instant, where Xei is closest to — c. Due to the ultrarelativistic behaviour, the current vanishes at these instants and 
we write: j{t) = ao^". After the calculations shown in appendix [A| the result can be expressed as 



£'cse(w) = f{uj) 



c(aia;)"+V2"+i 



(39) 



where gAi„(^) refers to a generalized Airy function, defined in Eq. (CI I and ^ 



Anyway, note that high order 7-spikes (n 3> 1) imply, that the nanobunch remains for a comparatively long time at 
low transverse currents. This appears to be unlikely here, as a static nanobunch would not stay together for long 
time without magnetic fields that can counteract the Coulomb explosion. Therefore, we go on to discuss only the two 
most likely special cases n = 1 and n = 2: 

1. The current changes sign at the stationary phase point and we can Taylor expand j(t) = aot. Consequently, 
Xei{t) = —vot + ait^ /3. The spectral envelope can now be written as: 



/(a;)(x|/>)|2c.-4/3 




(40) 



where Ai' is the Airy function derivative, ujrs ~ '^'^^^ y/o^lo, and 70 = (1 — fg) ^^"^ is the relativistic 7-factor 
of the electron bunch at the instant when the bunch moves towards the observer. As in the ROM models, the 



spectral envelope (40 1 does not depend on all details of the electron bunch motion Xei, but only on its behaviour 



close to the stationary points, i.e. the 7-spikes. 

2. In the case, when the current does not change sign at the stationary phase point, we Taylor expand j{t) = aot'^ 
and Xei{t) = —Vot + This yields to the spectral envelope 



/(w)(x|/>)pu;-6/5 




(41) 



with S" being the second derivative of S(a;) = gAi2(x) — {2n) ^ 

X / exp [i (^xt + t^/5)] dt, a special case of the canonical swallowtail integral For the characteristic fre- 

quency uJrs we now obtain oj^s ~ 2^/'*^/a77Q'^. Because now even the derivative of x^i is zero at the stationary 
phase point, the infiuence of acceleration on the spectrum decreases and the characteristic frequency scaling is 
closer to the 7^-scaling for a mirror moving with constant velocity. 



In Fig. [9] the CSE spectra of the synchrotron radiation from the electron sheets are depicted. Comparing them to 
the 8/3-power law from the BGP-case, we notice that, because of the smaller exponents of their power law part, the 
CSE spectra are much fiatter. E.g., around the 100th harmonic we win more than two orders of magnitude. Note that, 
as in the case of higher order 7-spikes in the ROM model, side maxima are found in the spectrum (41 ). This might 



provide an explanation for modulations that are occasionally observed in harmonics spectra, compare Ref. [27. .30j . 
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Figure 9: Coherent ID synchrotron spectra fo r an infinitely thin electron layer f{uo) = 1 and LUrs — 100. The blue line 
corresponds to Eq. (401 and the red line to Eq. (41 1. For comparison, the dashed black line denotes the BGP 8/3-power law. 




Figure 10: Normalized transverse fluid velocity Vy/c of the electron nanobunch. Parts of the plasma with a density below 
500 iVc are filtered out. 

To compare the analytically obtained spectrum with the PIC result, the finite size of the electron bunch must be 
taken into account. We assume a Gaussian density profile which leads us to 



|/H|2=CXp 



UJrf 



(42) 
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Figure 11: Dependence of the intensity boost r; — ma.x(E^)/ max{Ef) and the pulse compression F — {loot)~^ , where r is the 
FWHM width of the attosecond intensity peak in the reflected radiation, on ao- The laser amplitude ao is varied between 5 
and 195 in steps of 10. Other parameters are the same as in Fig. [Tja. 



Thus the spectral cut-ofF is determined either by w^s, corresponding to the relativistic 7-factor of the electrons, or 
by LUr f corresponding to the bunch width. A look at the motion of the electron nanobunch in the PIC simulation 
(Fig. 10 1 tells us that there is no change in sign of the transverse velocity at the stationary phase point, consequently 
we use Eq. (41 1. We choose ujrf — 225 and ujrs — 800 to fit the PIC spectrum, corresponding to a Gaussian 
electron bunch f{x) — exp [— (x/fJ)^] with a width of S = 10^'^A and an energy of 7 ^ 10. This matches well with the 
measured electron bunch width (5fwhm — 0.0015 A [see Fig. [7|^b)] and the laser amplitude ag — 60, since we expect 
7 to be smaller but in the same order of magnitude as ao. In this case Urf < ^rs, so the cut-off is dominated by 
the finite bunch width. Still, both values are in the same order of magnitude, so that the factor coming from the 
Swallowtail-function cannot be neglected and actually contributes to the shape of the cut-off. The modulations that 
appear in Fig. [9] for frequencies around uij-s and above cannot be seen in the spectra, because it is suppressed by the 
Gauss-function Eq. ( |42[ ). The analytical synchrotron spectrum agrees excellently with the PIC result, as the reader 
may verify in Fig. [6]jb) . 



3. Sensitivity of the nanobunchmg process to parametric changes 



Now, we have a look at the dependence of the harmonics radiation in and close to the nanobunching regime on the 
laser and plasma parameters. Exemplary, the laser intensity and the pre-plasma scale length are varied here. The pulse 
duration however will be left constantly short, so that we can simply focus our interest on the main optical cycle. For 
longer pulses, the extent of nanobunching may vary from one optical cycle to another, which makes a parametric study 
more difficult. We are going to examine two dimensionless key quantities: the intensity boost 77 = max(_E^)/ m.ax{Ef) 
and the pulse compression F = {ljqt)^^. It is straightforward to extract both magnitudes from the PIC data, and 
both are quite telling. The intensity boost 77 is a sign of the mechanism of harmonics generation. If the ARP boundary 
condition Eq. (14 1 is approximately valid, we must of course have 77 « 1. Then again, if the radiation is generated 



by nanobunches, we expect to see strongly pronounced attosecond peaks [see Eq. (49l[ in the reffected radiation 



and therefore 77 3> 1. The pulse compression F is defined as the inverse of the attosecond pulse duration. In the 
nanobunching regime, we expect it to be roughly proportional to rj, as the total efficiency of the attosecond pulse 
generation remains T^atto ^ 1, compare Eq. (47 1. In the BGP regime, there are no attosecond pulses observed without 
spectral filtering. So the FWHM of the intensity peak is on the order of a quarter laser period, and we expect F ^ 1. 

In figure [TT] the two parameters 77 and F are shown in dependence of oq. Except for the variation of oq, the 
parameters chosen are the same as in Figs. |6] [Tjb) and[lO| 



First of all we notice, that for all simulations in this series with ag ^ 1, we find 77 ^ 1. Thus, Eq. (14 1 is violated 
in all cases. Since we also notice F ^ 1 and F '-^ 77, we know, that the radiation is emitted in the shape of attosecond 
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Figure 12: Dependence of the intensity boost i) — max(ii'^)/ max(i5f ) and the pulse compression F — (cjor)~^, on the plasma 
scale length in units of the laser wavelength L/A in the lab frame. Except for the plasma scale length, parameters are the same 
as in Fig. [Tb. The plasma ramp is again an exponential one oc exp(a;/L). 



peaks with an efficiency of the order 1. This indicates, that we can describe the radiation as CSE. The perhaps most 
intriguing feature of Fig. (TTjis the strongly pronounced peak of both curves around = 55. We think that because of 
some very special phase matching between the turning point of the electron bunch and of the electromagnetic wave, 
the electron bunch experiences an unusually high compression at this parameter settings. This is the case that was 
introduced in subsection IIII D ll 

Figure [12] shows the two parameters rj and T as functions of the plasma gradient scale length L. It is seen that both 
functions possess several local maxima. Further, rj and F behave similar apart from one runaway value 'At L — 0.225A, 
where the FWHM peak duration is extremely short, but the intensity boost is not as high. A look at the actual field 
data tells us that in this case the pedestal of the attosecond peak is broader, consuming most of the energy. This 
deviation might e.g. be caused by a different, non-Gaussian shape of the electron nanobunch. 

The maximum of both functions lies around L = 0.33A, the parameter setting analyzed in detail before. In the 
limit of extremely small scale lengths L < O.IA, rj and F become smaller, but they remain clearly bigger than one. 
Thus the reflection in this parameter range can still not very well be described by the ARP boundary condition. For 
longer scale lengths L > G.8A, both key values approach 1, so the ARP boundary condition can be applied here. This 
is a possible explanation for why the BGP spectrum (27 1 could experimentally be measured at oblique incidence [20J. 



E. Harmonics emission in forward direction 



Up to now, we have discussed the harmonics emitted at the front side of an overdense foil, propagating in backward 
direction together with the reflected light. If however the foil used for HHG is sufficiently thin, harmonics are also 
emitted in the forward direction, albeit to a weaker extent |31H33j . 

Two main mechanisms [33] can be made responsible for the production of harmonics at the back side of the foil: 

1. Fast (Brunei) electrons which are transmitted through the foil may trigger the emission of harmonics not only 
at the front side of the foil, but also at its backside. Here, harmonics are produced up to the maximum plasma 
frequency as is characteristic for the sub-relativistic regime. This kind of radiation depends strongly on the 
properties of the density gradient on the backside. For a too steep gradient, only the very weak transition 
radiation is generated, whereas for the right scale length, harmonics are observed due to the CWE mechanism. 

2. Compressed electron bunches at the front side can emit high frequency synchrotron radiation in both directions. 
Radiation with frequencies above the maximum plasma frequency is transmitted through the foil, resulting in 
forward emission of high harmonic orders. We focus on this second mechanism here, as it is dominant in the 
relativistic regime. Obviously, the "oscillating mirror" model has no business here, since we are talking about 
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Figure 13; Electron density and contour lines (cyan) of the emitted high frequency radiation (a; > ujp). Relativistic harmonics 
emission is observed in both forward and backward direction, both generated at the front side. Simulation parameters are: 
laser amplitude ao — 10, 45° p-polarized incidence, maximum plasma density: Ue = 100 ric, exponential density profile at the 
front side with scale length L = 0.06 A and steep edge at the rear side. 



the transmitted, not the reflected part of the light. 

Let us also note, that the characteristics of the harmonics emitted in the forward direction may well be different 
from the ones emitted in backward direction, as two different 7-spikes are responsible for their generation. 
Whereas high frequencies in the backward direction are emitted, when the electrons are moving with maximum 
velocity away from the surface, they are emitted in forward direction when the electrons are moving with 
maximum velocity towards the surface. Thus the generation efficiency of forward harmonics depends on the 
existence of such a 7-spike, and on the compression of the electron bunch during this instance. 

An example of the second mechanism is shown in Fig. |13[ It is seen, that in this case both forward and backward 
harmonics are generated at the front side of the target. Further we observe (see detail in Fig. 13 1 that they are 
generated by two distinct electron bunches. The backward harmonics are mainly emitted by a bunch propagating 
away from the surface and becoming dispersed when returning into the plasma. In contrast, the forward harmonics 
are generated by a second bunch that achieves its highest compression when travelling into the plasma. 

Fig. [14] displays the spectra of the radiation emitted in forward and backward direction from the same simulation 
run. We notice, that the spectrum of the forward harmonics does not contain much radiation at low frequencies 
uj < LOp. This is easily understood due to the fact that the harmonics generated at the front side have to propagate 
through the plasma. Harmonics generation at the rear side due was suppressed in our simulation by the sharp density 
edge at this side. Further it is seen that the spectra of forward and backward emitted harmonics fall off at different 
rates. This is no surprise, since they are generated at different instants and even by different electron bunches, as 
seen in the detail graph inside Fig. |13[ At their respective 7-spikes, they possess distinct compression, shapes and 
energy, leading to different spectral slopes. 



IV. RELATIVISTIC HARMONICS RADIATION AS ATTOSECOND PULSE TRAIN 



A. Characterization of the generated attosecond pulses 



As the calculations above have shown, the emitted harmonics are phase locked. Therefore, they form a train of 
extremely short pulses. The shortest possible pulse duration can be estimated by T ~ l/w^, wherein is the 
characteristic roll-off frequency of the harmonics generation process. If the spectral roll-off is due to the relativistic 
7-f actor and the 7-spike is of the first order, this means that 



(43) 
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Figure 14: Spectra of radiation emitted in forward (red) and backward (cyan) direction, taken from the same simulation as 
Fig. [13] The spectrum of the incident laser (blue), the plasma frequency (black dotted) and an 8/3-power law (black dashed) 
are also charted for reference. 



Thus for the typical values of 7, the achievable pulse duration is in the order of a few attoseconds or even less. 
Note that the cubic scaling in 7 exceeds the possible pulse compression by the simple Doppler effect, which yields a 
duration not shorter than T^oppier ~ 1/(47^)- 

The 1/7^-scaling can physically be understood by having another look at a characteristic 7-spikes. For the ROM- 
process, this is 7ARp(i) as defined in Sec. IIIB3 and for the coherent synchrotron emission this is the 7-factor 
corresponding to the longitudinal bunch velocity component Xeiit). Assuming a 7-spike of first order here, we have 
v{t) w Do — awgi^ around the maximum. Consequently, the 7-factor can be written as 
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Evaluating the temporal width of the spike in ^{t) at t = yields Ai ^ 1/ {uiq^qO^/'^). Since the high order 
harmonics are produced only during the 7-spike, the duration of the corresponding attosecond pulses are in the same 
order of magnitude, in agreement with (43 1. 



1. Attosecond pulses from ROM harmonics 

In order to unravel the attosecond pulses contained in the ROM harmonics radiation, it is required to filter out the 
lower harmonic orders. The high-frequency cutoff of the power-law spectrum defines the shortest pulse duration that 
can be achieved this way. 

Assuming that the harmonics are emitted coherently and in phase, we expect their duration Tatto to be roughly 

the inverse of the absolute spectral width (ASW) . So as to estimate the pulse duration achieved by a certain filter 

1/2 

we therefore calculate the ASW Aco = ^(w^) — (w)^^ . This can be done analytically for a typical BGP spectrum. 
The spectral high-pass filter can be introduced as a sharp low- frequency cutoff beyond cjiow We write 



Iq lu exp(— w/wj.) ^('i^—'^iow), substituting the Airy function in Eq. (27l by an exponential for the sake of simplicity. 
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As shown in Ref. [S], this is a reasonable approximation here. Then we obtain for the ASW: 



(Ac) 



Ti^l; x)rC^; x) - [Ti-I; x)]' 



(45) 



wherein r(s; x) = t^^^ exp(— t) dt is the upper incomplete gamma-function and x = ujiow/^r- 
The function, normalized to ojr, is plotted in Fig. [15] For very low filter frequencies wiow ^ <^r, the ASW is still 
much smaller than Ur, consequently the observed attosecond pulses are not as short as they can be according to 
Eq. (43 1. If the filter frequency is increased up to ujr or even beyond, the ASW tends towards Alu — >■ Ur, therefore 
the generated pulses approach the duration given by Eq. ( 43 ) . 

Another interesting property can be extracted solely by looking at the spectrum. If the filter frequency is well 
below the critical roll-off frequency wiow ^ w^, the power law part of the spectrum dominates, resulting in a large 
relative spectral width (RSW) Alo/uj > 1. Consequently, single-cycle pulses are to be expected. If on the other 
hand wiow 2> w^, the spectrum decays exponentially above the filter frequency, leading to a relatively small RSW 
Auj/uj <C 1. Analytically we obtain for a BGP spectrum: 
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Again, the function is plotted in Fig. 15 



It is seen that Aw/w — >■ for a; — >■ cx), so for high filter frequencies the 
spectrum becomes increasingly monochromatic, yielding many-cycle pulses. For comparatively low filter frequencies 
X — we find the asymptotic approximation (Aw/w)^ ss 4r(|)/(15 x^/'^) — 1, leading to an increasing relative spectral 
width in agreement with our previous thoughts. 

These analytical results compare well to the numerical ones shown in Fig. |16[ For wiow = 9.5 ~ 0.1 a;^ the pulse 
duration is approximately Tatto ~ 0.05 A/c = 5Tr, where = 21: /wr = 0.01 A/c. This agrees reasonably with what 
can be expected from Eq. (45l: Aoj k, 0.15^;^ for x = 0.1. Further notice that the pulse is single-cycle, in agreement 
with the large RSW Aw/ (tJ) « 0.77, as given by Eq. (|46]). 

Things look differently for the higher filter frequency wiow = 49.5 Wq ~ 0.5 w^. Now the pulse duration is about 
Tatto ~ 0.02 A/c = 2r,., roughly a little more than twice as short as before. This is in agreement with Aw w 0.34 w,- 
for x = 0.5 from Eq. (45 1. We also see that the pulse now contains slightly more than one optical cycle, indicating 
the somewhat smaller RSW Aw/ (w) « 0.44 given by Eq. (46). 
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Figure 16: Attosecond pulses in ROM harmonics radiation can be revealed by spectral high-pass filtering, (a) shows the 
spectrum of the radiation in a logarithmic representation, along with the used filter frequencies as red lines, the black dashed 
line represent the filter frequencies, (b) and (c) display the resulting attosecond pulses for the filter frequencies tJiow = 9.5 uiq 
respectively toiovt = 49.5 ujo- Here, the black dashed line denote an envelope computed as the absolute square of the corresponding 
analytic signal. 



2. Attosecond pulses from electron nanobunches 



As we see in Fig. [oja), the CSE radiation is emitted in the form of a single attosecond pulse whose amplitude is 
significantly higher than that of the incident pulse. This pulse has a FWHM duration of 0.003 laser periods, i.e. 9 as 
for a laser wavelength of 800 nm. This is very different from emission of the ROM harmonics, which need to undergo 



diffraction (see also Sec. VI I or spectral filtering [6j before they take on the shape of attosecond pulses. 

When we apply a spectral filter in a frequency range (wiow, i^high) to a power-law harmonic spectrum with an 
exponent q, so that I{uj) = Iq{ujo/uj)'^, the energy efficiency of the resulting attosecond pulse generation process is 
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(47) 



The scaling (|47| gives r]f^^ - (wq/wiow)^/^ for the BGP sped rum with q = 8/3. For unfiltered CSE harmonics with 
the spectrum q — A/3 the efficiency is close to r/^^^ = 1. This means that almost the whole energy of the original 
optical cycle is concentrated in the attosecond pulse. Note that absorption is very small in the PIC simulations shown; 
it amounts to 5% in the run corresponding to Fig. [6] and is even less in the run corresponding to Fig. [T] 



25 



The ROM harmonics can be considered as a perturbation in the reflected signal as most of the pulse energy remains 
in the fundamental. On the contrary, the CSE harmonics consume most of the laser pulse energy. This is nicely 
seen in the spectral intensity of the reflected fundamental for the both cases [compare Figs.[l];b) and[6];b)]. As the 
absorption is negligible, the energy losses at the fundamental frequency can be explained solely by the energy transfer 
to high harmonics. We can roug hly estimate this effect by I^^^/I^^^ « oj-^^^ du/ uj^'^^^ duj = 5. This value 

is quite close to the one from the PIC simulations: I^^'^' ^V4^'^' = 3.7. 

Further, we can estimate amplitude of the CSE attosecond pulse analytically from the spectrum. Since the harmonic 
phases are locked, for an arbitrary power law spectrum I{uj) cx uj~'^ and a spectral filter (wiow, i^high) we integrate the 
amplitude spectrum and obtain: 
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(48) 



Apparently, when the harmonic spectrum is steep, i.e. q > 2, the radiation is dominated by the lower harmonics cjiow 
This is the case of the BGP spectrum q = 8/3. That is why one needs a spectral filter to extract the attosecond 
pulses here. The situation changes drastically for slowly decaying spectra with q < 2 like the CSE spectrum with 
q — A/3. In this case, the radiation is dominated by the high harmonics Whigh- Even without any spectral filtering 
the radiation takes on the shape of an attosecond pulse. As a rule of thumb formula for the attosecond peak field of 
the unfiltered CSE radiation we can write: 
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Using rric = uJc/luq = 225, the lower of the two cut-off harmonic numbers used for comparison with the PIC spectrum 
in Fig. l6|b), we obtain -Epoak = 8.8 i?o- This is in nice agreement with Fig. l6[a). 



B. Isolation of single attosecond pulses by polarization gating 

Many applications in the field of imaging and control of quantum dynamics on the attosecond timescale [31 [3S] 
require single attosecond pulse instead of a pulse train. The single pulse can in principle be produced using a phase- 
stabilized single cycle laser. However, relativistic harmonics require a laser pulse intensity / 3> 10^* W/cm^ and pulses 
in this intensity range usually are several cycles long, leading to the production of longer attosecond pulse trains. 
Therefore, we are in need of a method to isolate single attosecond pulses from the pulse train. 

It was shown above that the attosecond pulse are emitted when the tangential components of the surface electron 
momentum vanish. This property can be used to control the HHG and to gate a particular attosecond pulse out 
of the train, see also Ref. pQ. In the ID geometry, the transverse generalized momentum is conserved: p_L — 
eAj^/c + pj^,o, where pj^ and Aj^are the tangential components of the electron momentum p and the vector potential 
A. Consequently, the attosecond pulses are emitted when the vector potential is zero. If the vector potential vanishes 
at several moments, there are several 7-spikes and correspondingly, several short pulses are observed in the refiected 
radiation, see Fig. [l7[a). To select a single attosecond pulse, we must ensure that the vector potential A_l turns zero 
exactly once. Since Aj^has two components, how often it vanishes depends on its polarization. For linear polarization 
under normal incidence it vanishes twice per laser period, while for elliptic polarization it never equals zero. A laser 
pulse with time-dependent polarization can be prepared in such a way that its vector potential turns zero just once. 
A pulse of time-dependent polarization can be equivalently represented as a superposition of two perpendicularly 
polarized pulses, driving and controlling pulse, with slightly different frequencies and phases. Our PIC simulations 
suggest that a controlling signal with a small fraction of the driver intensity is sufficient to manage the HHG, if the 
phase difference between the two laser pulses is chosen carefully. 

An example of this is shown Fig. |17[ b). In addition to the main pulse, which is the same as in panel (a) of the 
figure, a smaller controlling signal was used at an amplitude of oi = 6. Frequency and phase of the controlling pulse 
were slightly detuned in comparison to the main pulse. It is clearly seen that all attosecond peaks except for the 
main pulse at ct « 9.3 A are strongly attenuated. Only two distinctly smaller side peaks are left at ct « 8.8 A and 
ct « 9.8 « A. 

We conclude that it is possible to isolate single pulses in relativistic HHG via a polarization gating technique, even 
if the driver pulse contains more than one optical cycle. 
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Figure 17: Attosecond pulse train (a) without and (b) with the use of polarization gating technique. Parameters are: ao = 20, 
plasma density Ne — 90 Nc with a sharply defined surface. In (b) there is a second pulse with amplitude ai — 6 and polarization 
direction orthogonal to the first one, detuned to uji = 1.25 a;o and dephased by A(j> ~ it/8. 



V. LINE STRUCTURE IN RELATIVISTIC HARMONICS SPECTRA 



In section III theoretical models of surface HHG were discussed. From these models we were able to compute the 
envelope of the harmonic spectrum, but they do not tell anything about the structure of the individual harmonic 
lines. The line structure provides additional details about the laser-plasma interaction on the femtosecond timescale 
and thus may serve as a useful diagnostic. However, to utilize it, a thorough understanding is needed at first. This 
section aims to provide this understanding. 

In subsection |V A| we briefly discuss the line structure occurring in the moderately relativistic regime. In the 
highly relativistic regime, the spectral line structure is closely related to the phase of individual attosecond pulses 
inside the generated pulse train. Therefore, we examine the dependence of this phase on laser amplitude and plasma 
density in subsection |V B Next, we relate this to the chirp of the relativistic harmonics (Sub. VCl and calculate its 
spectral footprint (Sub. |VD ), which is well accessible in experiments. Such experiments have been conducted at the 
ARCTURUS facility in Dusseldorf. In subsection |VE[ we report about how they substantiate the presented theory. 



A. Spectral line structure in the moderately relativistic regime 

At moderate intensities, modulations in the spectral line structure such as half integer harmonics are mainly caused 
by parametric instabilities in the underdense part of the pre-plasma. 

Parametric instabilities, such as stimulated Raman scattering and the two plasmon decay in the underdense pre- 
plasma lead to creation of plasmons at about half the laser frequency |23| . These plasmons can then recombine 
with the laser or harmonics photons via sum frequency mixing, leading to side bands or spectral lines at half-integer 
multiples of the fundamental |36| I37| . This mechanism is prevalent for moderate intensities oq ~ 1, longer pulse 
durations cr ^ A and extended pre-plasmas. 

Moderate broadening of the harmonic lines may also be caused by the inherent chirp of the CWE process, see 
Ref. inni- This chirp arises due to the dependence of the excursion times of the Brunei electrons. For higher intensities, 
the excursion times are longer, thus the attosecond pulses are emitted with a longer delay. Assuming a bell shaped 
temporal profile of the laser pulse, this leads to a negative (blue to red) chirp. 
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At higher intensities oq ^ 1, the relativistic ponderomotive force of the laser sweeps away all electrons from 
the underdense plasma regions. Therefore, parametric instabilities play no important role anymore. Also, the CWE 
mechanism looses importance as the relativistic effects take over. However, for these pulses, there is again a mechanism 
that leads to a variation of the phase of the attosecond pulses depending on the temporal variation of the laser 
intensity. This can lead to heavy broadening and modulation of the harmonic lines, particularly for extremely short 
pulses CT > A. Let us now go on to discuss this mechanism in detail. We begin by numerically computing the 
dependence of the phase of the attosecond peaks on the laser intensity and other parameters. 



B. Attosecond peak phase in the highly relativistic regime 



This subsection is divided into the investigation of normal incidence and the investigation of s- and p-polarized 
oblique incidence. 



1. Universal phase relation in normal incidence 



We start by examining the case of normal incidence on a perfectly steep plasma boundary. To begin with, a suitable 
definition of the "phase of the attosecond pulse" is needed. 

Having another look at Fig.[lja), showing a quite typical case of the reflected electric field in normal incidence HHG, 
tells us what to do. Due to the discontinuities in the function Er{t), the time derivative possesses clearly pronounced 
peaks. Therefore, we define the "attosecond phase" <j) as the position of the maximum of the time derivative of the 
reflected electric field dfEr- Later on (Sec. VI I we will see, that dtEr also happens to play an important role in the 
computation of the far field. is normalized in a way, that = if there is only the Guoy phase shift in the case 
of simple non-relativistic refiection from an infinitely dense surface. With this definition, we measured for a huge 
range of densities = 20 . . . 450 A'c and laser amplitudes gq = . . . 450. The result is displayed in Fig. 18 
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Figure 18: Phase dependence of the attosecond peak on the inverse ^-parameter, S"^ = aoNc/N, under normal incidence. 

Very short pulses (r — 2tt/uj) were used, so that the surface remained intact during the interaction even for high 
intensities. The phase of the incoming laser pulse is chosen in a way so that Ei = Q a,t the maximum of the envelope, 
thus the attosecond peak is located close to the maximum of the envelope. The attosecond phase is plotted against 
the inverse S'-parameter (see Ref. |26i) S^^ ~ a^Nc/N . For all simulations in the highly relativistic regime oq ^ 1, 
we find an excellent agreement with the fit 




= 2.7 S'^i - 0.32, (50) 

while in the low intensity limit S^^ — > the phase shift tends to the value = acot((A^/A'c — 2) / {2^J N /Nc ~ 2)), 
which is expected from non-relativistic optics, approving the correctness of the PIC calculations once again. 

Physically, the phase shift </> can be understood as a consequence of the electron surface being pushed inside the 
plasma by the laser. If the electron surface is pushed in to a depth of A, we expect the phase to experience an 
additional shift cx A. Let us devise a rough model in order to understand the linear scaling of with . Therefore 
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Figure 19: Phase dependence of the attosecond peak on the inverse S'-parameter in the simulation frame, S^^ = aoNc^'' /N'-^\ 
under s-polarized oblique incidence; density A'^ and laser amplitude ao are varied. In relation to the laboratory frame S- 
parameter, Ss scales as Ss ~ Sl/ cos^ 6. Both the phase of s-polarized and p-polarized generated harmonics is displayed. The 
angle of incidence is 6^ = 45° 

we assume that there is a pressure balance between the ponderomotive force fpond Uq of the laser and the electrostatic 
restoring force fstat — qE cx TV^A^ of the plasma. Equalizing both terms yields A ^ gq/N and consequently, a linear 
dependence of </> on S^^. 

Note further, that in the ultrarelativistic regime the function (j) is indeed completely independent of the absolute 
plasma density. This is the clearest footprint of the S'-similarity |26j in laser-overdense plasma interaction observed 
so far. 



2. Phase behaviour at oblique incidence 

When considering oblique incidence, the polarization is crucial. For s-polarized oblique incidence, we retain a 
behaviour similar to the one observed under normal incidence. For p-polarized incidence, the behaviour changes in 
many ways. We analyze both cases using ID PIC simulations in a Lorentz transformed frame (see App. [B|. 

Consider Fig. |19[ As in the case of normal incidence, we confirm the dependence on the 5-parameter with high 
accuracy in the ultrarelativistic regime. If the density is varied, but the ratio S^^ = a^Nc/N is kept constant, there 
is no change in the attosecond phase. We also see that there is virtually no difference between the phase of the 
p-polarized and the s-polarized generated harmonics. This is evidence that they both are generated due to the same 
physical mechanism. They are not generated at separate phases as are CWE and ROM pulses in the weakly relativistic 
regime fT9| . 

Further, as in the normal incidence case, an approximately linear dependence on is found. This can be 
understood, as the mechanism leading to the indention of the electron plasma surface is basically the same as for 
normal incidence: There is a pressure balance between the ponderomotive light pressure and the electrostatic force. 

In the laboratory frame however, the ponderomotive light pressure is expected to be weaker compared to normal 
incidence, since the laser does not hit the surface head on, but under an angle 0. Seen in the simulation frame, the 
ions and the electrons possess currents in opposite directions. This generates a magnetic repulsion, counteracting the 
electrostatic restoring force. Effectively, it leads to a mitigation of the electrostatic force by a factor of I/7. Therefore 
we expect, that the scaling in s-polarized oblique incidence should remain independent of the angle 9 if we consider 
it a function of ^eff = Ss/j — Sl/ cos^ 6. 

This can well be confirmed by the numerical results depicted in Fig.|20] We conclude, that for s-polarized incidence, 
the phase of both the s-polarized and the p-polarized fraction of the generated harmonics is determined only by the 
effective S'-parameter SeS = aoNc/{N cos^ 9) and does not depend on ag, N and 9 separately. 

For p-polarized incidence, matters are more complex. In addition to the ponderomotive force, the surface is also 
pushed in and pulled out directly by the longitudinal electric field component of the laser. Numerical results are 
shown in Fig. [2l] 

Despite of the highly complex interaction, the attosecond phase (f) again depends only on the S'-parameter, not on 
ao and N separately. The slope is however not linear anymore. As we can see from Fig. |21[ the non-linearity increases 
with the angle of incidence 9. 
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Figure 20: Phase dependence of the attosecond peak on the inverse effective S'-parameter under s-polarized incidence, angle of 
incidence 6 and laser amplitude ao varied. The effective 5'-parameter is defined as Sefi = Sl/ cos^ 9 = Ss/y to bring out the 
consistent linear dependence. Here, only the phase of the s-polarized generated harmonics is displayed since the one for the 
p-polarized harmonics almost agrees (see Fig. 19 1. 
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Figure 21: Phase dependence of the attosecond peak on the inverse S-parameter in the simulation frame under p-polarized 
incidence, angle of incidence 9, density A*' and laser amplitude ao are varied. 



Note further that in the case of p-polarized incidence, the duration of the pulse may also play an important role. 
Oblique p-polarized incidence can lead to the generation of very strong quasi-static magnetic fields close to the surface. 
Therefore, memory effects are present and (/) is not just a function of the instantaneous intensity but a functional of 
the whole history of the incident field. In section [V E| we will show an example of this highly interesting effect. 



C. Evidence of harmonic chirp in PIC simulation 

As we have just seen, the phase of the attosecond pulses generated from overdense plasmas depends on the S*- 
parameter of the interaction. Because of the relativistic radiation pressure of the laser pulse, the electrons are pushed 
inside the plasma during the rising edge of the laser pulse, causing an initial red shift of the reflected light. Later, the 
electron fluid will return to its original position and therefore cause a blue shift. This shifting of frequencies is called 
harmonic chirp[49_. 

The motion of the reflecting surface can be followed in an x-i-colourscale image of the transverse magnetic field 
component. Figure [22] shows such an image for a PIC simulation of normal laser incidence on a perfectly sharp 
plasma boundary. Realistic cases with oblique, p-polarized incidence on a plasma with a finite density gradient will 
be discussed in subsection |VE[ In our simple case we see that the motion of the surface is well described by a Gaussian 
function, i.e. the surface displacement is proportional to the instantaneous laser amplitude. This agrees with the 
observations made in subsection |V B| with even shorter laser pulses. 
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Figure 22: Transverse magnetic field component of a laser reflecting and generating harmonics at an overdense plasma 
surface at normal incidence. The dashed line highlights the motion of the reflecting surface, fitted by a Gaussian function, 
corresponding to Eq. ( |56[ ) with a — 1.35. Parameters are: laser amplitude ao = 10, duration cr — 5A; sharp edged plasma with 
density rie = 20nc starting at a; = 1 A, fixed ions. 
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Figure 23: Spectrogram with logarithmic colourscale of the refiected radiation. Simulation parameters are the same as in 
figure |22| For the computation of the spectrogram, a Blackman-Harris window with a width of about eight laser periods was 
used. 



As a result, the reflected radiation contains a positive chirp. This chirp can be made visible in a time-frequency 
image (or spectrogram). To compute the spectrogram, the time-series data from the PIC simulation is multiplied with 
a bell-shaped window function that is gradually moved over the data. Then, spectra of the products are calculated, 
yielding the spectrogram. 

In Fig. |23] such a spectrogram is shown. The data stems from the same simulation as the previous figure. The chirp 
from red to blue is clearly visible in the fundamental and all harmonics shown. The higher the harmonic number, the 
more pronounced is the chirp. This is due to the nature of the Doppler effect. Because the relative Doppler frequency 
shift Alo/lo is constant, the absolute shift Auj is proportional to the frequency. Further note that the positive chirp 
(from red to blue) observed here distinguishes this relativistic regime from the CWE regime, where a negative chirp 
(from blue to red) is observed ^38) . 
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D. Spectral footprint of harmonic chirp 



In the time integrated spectrum the chirp is visible in the Hne structure. Due to the dependence of the Doppler 
shift on Lo, we also expect the line structure to vary according to the spectral region. 

To describe this structure analytically, assume the radiation is given as a sequence of attosecond pulses emitted at 
the times tk with identical shape f{t) but possibly different amplitudes Ek, corresponding to the laser envelope: 

= ^i;fe/(i-ifc). (51) 



Now we Fourier transform Eq. (51 1 and take the absolute square to arrive at the spectrum: 



(52) 



where /(w) denotes the Fourier transformation of the attosecond pulse shape function f{t). Its absolute square |/(aj)p 
corresponds to the spectral envelope that has been discussed in Sec. |III| Here, we concentrate on the second factor 
J(a;) that represents the spectral line structure. 

In the trivial case of equidistant pulses with constant intensities, i.e. E^ = 1 and tk = kT^ = iirk/LOQ, the result is a 
sequence of sharp harmonic lines at multiples of the fundamental frequency . Such a spectrum occurs for harmonics 
generated by comparatively long laser pulses (picosecond range) with moderate intensities. Early experiments on 
surface harmonics generation worked with such pulses and obtained spectra close to this prediction, compare e.g. 
Ref. [81. 

The spectrum changes as pulses become shorter and more intense. As described above, the harmonics move from a 
red-shifted to a blue-shifted phase due to the Doppler effect of the averaged surface motion. To get a first impression 
of the effect on the spectra, let us consider two trains of pulses with a slightly different periodicity Ti and T2 . Both 
pulse trains will produce a train of harmonic lines corresponding to their repetition frequencies uii = 2n /Ti (i G {1, 2}). 
The harmonic lines will then interfere with each other. Provided they possess a finite linewidth Suj and the difference 
between the two interfering frequencies is small in the sense A = a;2 — wi <C we can calculate the frequency 
period fl of the occurring interference pattern by setting fl = nuj2 = {n + l)wi and therefore « Wq/A, where 
Wo = (wi + u)2)/'2 is the centre frequency. 



Fig. 24 shows the actual spectrum J(w). The figure confirms the analytically evaluated modulation frequency of 
= 10 loq. The first few lines clearly pronounce multiples of the fundamental frequencies, then the lines become 
broader and at around Sajg, we observe two lines per harmonic. Around 10 wq, we observe one peak per harmonic 
again, but this time at half integer frequencies. After that, the structure repeats, shifting back to integer harmonics 
around 20 wq. 

Certainly, in reality the period of the as-pulses does not change abruptly, but continuously. Let us therefore consider 
a linearly chirped train of attosecond pulses: 

tk = '^-^(k-^kA, (53) 
where k = — m . . .m. Note that the parameter /3 represents the maximum cycle averaged velocity acquired by the 



reflecting surface. Inserting (53 1 into (52 1 we arrive at 



> Ek exp(-iujtk) = 1 + 2 > cos exp i . (54) 

The resulting spectra are depicted in Fig. |25] Again, we see a periodic modulation of the spectrum. Because of the 
continuously changing period lines also broaden, especially in the high frequency range. This leads to a decay of the 
modulation amplitude. Further, a quasi-continuum is observed at the spectral regions where two harmonics per period 
were observed with the discrete frequency model. The modulation period corresponds to the frequency difference A 
between the extreme ends of the linear chirp and therefore is proportional to the maximum surface velocity. Thus, it 
is possible to extract physical information from the line structure. If the reflecting surface attains a velocity of /3 (in 
units of c), we can expect large scale modulations in the line structure with a period Q given by 
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Figure 24: Spectral line structure of as-pulse train with two different frequencies lvi = 0.95 ojo and LO2 = 1.05 cjq- The pulse 
train contains 9 pulses. 
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Figure 25: The model from Eqs. ([52| and 1 53 1, using 1/(1^)1 = 1 and m = 64 throughout. Spectra are smoothed by convolution 
with a Gaussian function of FWHM Auj — 0.1 a;o- 



Let us now reconsider the example from the previous subsection |V C[ From Sec. |VB| and also Fig. |22]of this section, 
we learn that the phase of the reflected radiation depends roughly linear on the laser field. Therefore, we can directly 
relate Ek and tk to the envelope g{t) of the laser. At normal incidence we further expect two attosecond pulses per 
period with alternating sign. This leads us to the model: 



tk ~ irk/cuQ + ag{'Kk/iUQ) 
Ek « i-lf g{TTk/iUo). 



(56) 



Inserting Eqs. (56l into Eq. (52) yields the structure of the spectral lines. For S = n^/aQric > 1, the model 



parameter a can be determined from the linear slope in Fig. 18 (or Eq. (50)). Therefore, a — 1.35 corresponds to 
S" = 2 and a = 0.54 to S' = 5. 
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Figu re 2 6: Comparison of the line structure from the fifth to the ninth harmonic (a) from PIC dat a an d ( b) due to Eqs. (52 1 
and (56 1. Simulation parameters corresponding to the blue line in (a) are the same as in Figs. 22 and |23[ the red line in (a) 
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Figure 27: Plasma density profiles corresponding to experiments at the Arcturus laser facility using AR and BK7 type plasma 
mirrors, as used in the PIC simulation. All magnitudes are shown in the Lorentz transformed simulation frame (see App. [B|. 



Now we can compare the simple model to spectra obtained from PIC data. As shown in Fig. |26] this comparison 
shows good qualitative agreement. Remaining differences can arguably be attributed to the non-linear dependence of 
the attosecond phase and the harmonics intensity in the only moderately relativistic interaction at the edges of the 
pulse, which are not included in the simple model. 

We conclude, that the line broadening observed in relativistic harmonics spectra can to a large extent be explained 
by the chirp due to unequal spacing of the attosecond pulses. It does not imply a loss in coherency of the individual 
attosecond pulses. 



E. Experimental confirmation of harmonic chirp 



Let us now have a look at a set of parameters derived from real experiments carried out at the Diisseldorf 
ARCTURUS laser facility. The laser is obliquely incident under an angle of 45° and the light is p-polarized with 
an estimated peak amplitude of about oq = 8. Two different kinds of plasma mirrors were used to improve the laser 
contrast ratio: the AR (anti-reflex coated, 0.1% reflectivity) mirror leading to a high contrast and an extremely steep 
density gradient and the polished BK7 glass (~4% reflectivity) yielding a medium contrast and a little less steep 
density gradient. In the PIC simulation, "double exponential" density profiles of the type 

, , J exp [a(x - xo)] -I- exp [6(a; - xi)] {x < X2) 
ne[x) = < ' (57) 
no {x>X2) 
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Figure 28: Colourscale image of the transverse magnetic field component as a function of time t and space x, overview and 
detail. The dashed lines in the detail graphs denote the part of the surface motion used to model the spectral chirp. The data 
stems from simulations with a realistic set of parameters: The incident laser has a Gaussian temporal profile a = ao exp{—t^/T'^ 



with amplitude ao = 8.1 and pulse duration r — lOA/c for both subfigures. The density profiles used are displayed in Fig. 27 
here (a) corresponds to the AR and (b) to the BK7 plasma mirror. 



were chosen as shown in Fig. [27] These profiles closely resemble the ones in the experiment, which could be estimated 
by simulations of the hydrocode Multi-FS[5 , conducted by Michael Behmke and Jens Osterholz. 

Figure [28] shows the transverse magnetic field from the interaction. Compare this to Fig. [22] The strong temporal 
asymmetry is conspicuous: Instead of instantly returning to its initial position as in the normal incidence case (Fig.)22|, 
the electron surface remains indented. The static magnetic field, created by the current of the Brunei electrons, holds 
the electrons inside. 

Because of this temporal asymmetry and the non-linear dependence of the phase on the amplitude, Eq. (56) ceases 
to apply here. To reproduce the exact spectral shape in our model would therefore require to exactly trace the surface 
motion with a non-elementary function. The purpose of our simple model is however not to exactly reproduce the 
spectrum, but to extract some crucial features. Our aim is to provide clear evidence that the modulations in the 
experimental and PIC spectra are caused by the unequal spacing between the attosecond peaks and to show, which 
information can be gained from the spectra. We therefore design the model as plain as possible, leaving only two free 
parameters /3 and to. We concentrate on the main phase of harmonic generation t = 30. . .45A/c and approximate 
the surface motion during this phase by a parabola, corresponding to a linear chirp. Also, the temporal asymmetry 
is ignored, taking the sum in Eq. (54 1 always from —to to to. 

Considering the complexity of the actual, highly non-linear process and the simplicity of the model, we observe an 
excellent agreement between the experimental spectra, the PIC spectra and the analytically calculated ones, compare 
Fig. [29] From the conjunction between the plasma motion in Fig. [28] and the model spectra from Fig. [29] it becomes 
clear that the modulations in the spectrum are caused by the unequal spacing between the attosecond pulses within 
the generated pulse train. 

Another interesting detail is the slight redshift to be observed in experimental [Fig. 29 b)] and PIC [Fig. 29 d)] 
data, in particular between the 15th and the 25th harmonic. This overall redshift is a footprint of the aforementioned 
temporal asymmetry in the femtosecond plasma dynamic. Therefore, it witnesses the self-generated static magnetic 
field. 

Let us now estimate the surface velocity from the experimental spectra, employing Eq. ( 55 1 . In the "medium 
contrast" case, we observe a transition from integer harmonics in the region up to the 25th order to half-integer 
harmonics in the region beyond the 30th order. In between, the lines merge into a quasi-continuum. Thus, the 
modulation cycle is about il ^ 27a;o, and the cycle averaged surface velocity is in the order of 0.01 c. 

We conclude that the harmonic spectrum is rich in information about the plasma dynamics on the femtosecond 
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Figure 29: Comparison of (a)-(b) experimental, (c)-(d) PIC and (e)-(f) model spectra according to Eqs. (52 1 and (53 1. The 
left column represents the experiment conducted with the AR plasma mirror (high contrast, very short pre-plasma) and the 
corresponding simulations, the right column represents an experiment conducted with the BK7 plasma mirror (medium contrast, 
slightly longer pre-plasma) and the corresponding simulations. Parameters in (c) and (d) are the same as in Fig. 28 Parameters 
chosen in model: (e) /3 = 0.0028, (f) (3 = 0.011, \fiuj)\ = o;"*/^ and m = 7 for both. PIC and model spectra were smoothed 
by convolution with a Gaussian function of FWHM Alj — 0.15 cjq for the sake of better comparability with the spectrometer 
data. Intensity units within one row are comparable, but non in between the rows. 



timescale. In the presented experiment, the influence of the laser contrast on the pre-plasma scale is clearly reflected 
in the harmonics chirp and thus, in the spectral line structure. From modulations in the spectrum, we can estimate 
the cycle averaged velocity of the electron plasma surface during its interaction with the main pulse. 



VI. RELATIVISTIC HARMONICS IN 3D GEOMETRY 



Up to now, we studied the theory of surface HHG in a simplified ID geometry. This chapter is dedicated to the 
investigation of harmonics spectra and the attosecond pulses in realistic 3D geometry. There are two new groups of 
phenomena to be considered: 

1. Due to the extremely broad spectrum of the generated radiation, diffraction will exert a major influence on the 
spectrum in the far field. Although the photon fiux through a surface surrounding the whole solid angle of 2tt 
in front of the plasma remains constant for each frequency, we expect the radiation field to be inhomogeneous 
so that the spectrum changes as a function of position. Firstly, because of the different diffraction lengths of the 
harmonics, and secondly, because the harmonics field distribution at the plasma surface may differ in intensity 
and phase from that of the driving laser. It is our main concern to investigate these phenomena since they may 
lead to very useful applications. 

2. In particular for very tightly focused laser pulses (cr ^ A), 3D effects may play a role in the physics of the 
harmonics generation itself, so that ID theory is not applicable anymore. We discuss the transition to this 
regime in section [VI D| 

Two works [391^0] precede the study of this section, which is based on the author's own work Ref. [¥T1H21. The letter 
[39] by Gordienko et al. studies the coherent focusing of the surface harmonics radiation under strongly idealized 
conditions, and Naumova et al. |40j concentrate on extremely tightly focused ("A'^-regime") pulses. Here, we present 
a broader overview of 3D phenomena that play a role in HHG experiments. 
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We begin with some analytical estimations about how the harmonics spectrum changes due to diffraction in vacuum 
(subsection VIA). Then (subsection VI B I. we verify the assumption of "locally independent" HHG, to see down to 
which focal spot sizes ID theory is still applicable for the generation process itself. We propose to exploit diffraction 
effects as "spatial spectral filters" (subsection |VI C I - an alternative or supplement to spectral transmission filtering to 
unveil the attosecond pulses in the harmonics radiation. The coherent focusing of harmonics (CHF), as proposed in 
Ref. j39] , is a promising pathway towards extremely strong fields, perhaps reaching up to intensities that should allow us 
to observe exotic effects such as non-linear vacuum polarization |43II44| . predicted by quantum electrodynamics (QED). 
In subsection |VI D| we discuss it, thereby considering a more realistic setup compared to the original proposition in 
Ref. 1311 • 



A. Harmonics spectrum changes due to vacuum propagation 

Vacuum propagation exerts an intriguing influence on high harmonics radiation generated from solid surfaces. We 
start with simple analytical estimations in order to illustrate this. Further, we aim to give an idea of how these effects 
might be harnessed to work in our favour. 

We begin by considering a linearly polarized Gaussian laser pulse normally incident onto a planar and sharp-edged 
overdense plasma surface. The generated harmonic spectrum can be approximated by a power law everywhere close 
to the surface, so that /near('', ^) ~ ^o('') (w/wo)^^ for lo ^ wq, where the exponent p depends on the exact HHG 



mechanism, e.g. p = 8/3 in the BGP case (see section IIIB3). If we neglect the intensity dependence of absorption 



/o(r) is proportional to the intensity of the incoming beam, and therefore Gaussian. 

Let us at first regard the idealized case that the phase of the generated harmonics does not depend on r. Since 
the frequency w is of course very different for distinct harmonic orders but the focal spot size a is the same for all of 
them, we easily find that the beam divergence = Oquio/oj is inversely proportional to the harmonic order. The high 
orders are emitted into a much smaller cone than the lower orders. On the optical axis, the spectrum in the far-field is 
therefore expected to show a much stronger pronunciation of high harmonic orders compared to the spectra obtained 
within ID models. 

Let us therefore evaluate the development of the spectrum on the optical axis. We are also interested in including 
the intensity-dependent frequency cut-off in our estimate, so we start with: 

/near(r-,c^) =/o(r) (^^^ " ([a(r)]« l., - w) . (58) 

Here, the cut-off frequency is assumed to have a power law dependence on the vector potential amplitude a{r) at a 
given point of the surface. 

The general starting point for our computations of the far-field is the Kirchhoff formula: 



in 



dA ■ (G(f, r')Vij{T^) - V^(f)VG(f, f')) 



(59) 



with G(f*, f) — exp(iu;|r — r'\)/\r — r*]. Specialized to the geometry of a beam focused onto a planar surface, assuming 
cylindrical symmetry and using F- = 0.5{Ey — B^) this becomes 



r dr dtF^ x = 0,r,t 



(60) 



where \x\ 3> Rmc 

calculations may easily be re-interpreted for curved surfaces, as we are going to see in section |VID 



a is assumed. Equation (60) is specialized to planar surfaces, but the results of the following 



Applying (60) to the harmonics generated by a Gaussian laser pulse as given by Eq. (58) and no phase dependence 



on r included, we find the far field spectrum to be 

I{x,Ll>) = Iq 



(<^o) 



flo V 



X. ao-^oo X 



(61) 
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wherein xri = ttct^/A is the Rayleigh length of the fundamental. 



Eq. (61 1 shows explicitly how vacuum propagation influences the harmonics spectrum on-axis. Just by picking the 
right point in space in front of the harmonics-generating surface, we can find a spectrum decaying two powers slower 
than the spectrum predicted by ID theory, i.e. (using the BGP exponent p = 8/3) / oc w^^/^^. Physically, the reason 
for this is the much stronger coUimation of the higher harmonic orders. 

Another interesting detail is that the sharp spectral cut-off in the near-field yields a soft roll-off in the far-field. 
Provided the far field spectrum can be measured accurately, conclusions on the constants q and ujc, which determine 
the general intensity dependence of the harmonics cut-off (see Eq. ( |58[ )) in the near-field are possible in principle. 

However, the above calculation presumed that there is no phase dependence on the distance from the optical axis 
r in the near field. In section |V B[ we have seen that the attosecond phase (p depends on the S'-parameter of the 
interaction. Therefore, to produce the effect explained above, it is necessary to keep the local S'-parameter constant 
along the surface. In section [VI C| we discuss methods to achieve this, employing PIC simulation to substantiate our 
proposal. 

Now let us estimate the consequences of the variation of the local S'-parameter along the surface in the case of a 
Gaussian laser pulse, normallv|50| incident on a surface with steep density gradient up to a constant density. Due 



to the curved phase surface, we expect the harmonics pulse to self-focus. Applying Eq. (50 1, we can calculate the 
curvature of the generated phase surface and consequently, the self-focusing distance x^f'- 

Xsf = ^^Ri- (62) 
Due to this self-focusing, the divergence angle of the individual harmonics is not simply proportional to the harmonic 



wavelength anymore. Equation (50 1 allows us to derive an expression for the divergence angle: 

(63) 




7A — Co , , ■ , r, 

wherein Oq is the solid angle, from which the laser itself is focused. 

B. Checking the assumption of locally independent generation 

In the above calculations we have assumed, that the generation process itself can be described by the ID models 
discussed in section [Till Only to investigate the diffraction of the emerging radiation, we consider the real 3D geometry. 
In other words, we presumed that the harmonics are generated locally independently at each point of the surface. This 
means that the radial field gradient has no infiuence on the harmonics spectrum and phase at a certain point. There 
is no transverse energy transfer. Mathematically this condition can be expressed as 

F_(y,z,<)-i^_[F+|^,^^ (64) 

where F_ stands for the refiected field and F+ for the field of the incoming radiation. 

Note that this assumption also allows us to perform ID instead of 3D PIC simulations, even if we are interested in 



the far-field of a realistic 3D geometry. Assuming the validity of Eq. ( 64 1 , we can merge the results of a series of ID 



simulations, utilizing Eq. (60 1 to obtain the far-field. ID simulations are computationally much cheaper and so they 



can be performed with a higher resolution in the same amount of time on the same computer. 



Let us now check in which parameter region the condition (64 1 is satisfied. The refiected field is of course generated 
by plasma electrons. The electrons are driven by the electromagnetic field of the laser pulse. In the ultra-relativistic 
case, the size of the electron orbits is on the order of A. The scale length, on which the radiation intensity at the 
surface changes in radial direction is the beam waist a. Therefore, if cr ^ A, the electrons might mediate between 



regions of different intensities, endangering the validity of Eq. (64 1. If cr ^ A, the electrons are not able to travel this 



distance and one expects that (64 1 is fulfilled. 

We now compare results of ID and 3D PIC simulations to verify this. 3D simulations were performed using the spot 
sizes cr = 5, 2, 1, 0.5 A, all with the dimensionless laser amplitude ag — 30 and the plasma density — 90 Uc, where 
Uc = ujQin/A'iTe'^ is the critical density. The laser pulse is linearly polarized in y-direction. The refiected field in the 
3-dimensional PIC simulations was always recorded at a distance of lA to the originally sharp-edged plasma surface. 
This field is directly compared to the result of ID PIC simulations with the same parameters {6id,3d)- Anyway, for 
the most tightly focused pulses the field recording distance is in the order of the fundamental Rayleigh length so that 
a direct comparison with ID PIC is pointless. Instead of this, it can be compared to the far field calculated with 
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Table II: Deviations from the assumption of locally independent harmonic generation according to the measure (651. 



equation (60 1 from ID results. For a = 1\, the distance is already too big to compare directly but yet too small to 
use (60 1, so the value is missing here. 

Yet two more comparisons were performed. 



The absolute far field (x — ?> oo) calculated by Eq. (60 1 from ID and 
This method verifies the accuracy of our quasi-lD calculations of the far-field 



3D PIC data is compared ((^i'^^sd)- 
directly. 

The third comparison {6y^z) concerns the radial symmetry of the reflected pulse. In the 3D geometry, it can in 
principle be broken because of the linear polarization of the incoming laser pulse, but it obviously cannot be broken 
Thus this symmetry check provides another indirect criterion for verifying (64 1. The fields 
An 2, z = 0, t) and F^{x ^ lA, y = 0, z = 



as long as (64 1 holds 



both obtained from the 3D 



compared are F^{x — lA, y = ay 
PIC simulations. 

The results of all these comparisons are collected in Tab. [IT] The relative deviations are measured using 



h,2 



!dtm)-Mtr 

fdt\h{tW + \Mtw 



(65) 



This measure is if /land /2 are identical functions and 1 if they are completely uncorrelated. 
It can be seen that deviations are very small for not too tiny focal spots such as cr = 5A. Thus, models based on 
Eq. (64 1 may be used for the vast majority of today's HHG experiments. Even in the case a = X, Eq. (64) still holds 



as a rough approximation. Our studies of even tinier focal spots have shown, that the deviations from Eq. (64 1 are 
generally not favourable for the generation of attosecond pulses. 

We now go on to discuss diffraction effects, assuming the validity of Eq. (64 1, which is correct for not too small 
focal spots in the sense explained above. 



C. Self-focusing and spatial spectral filtering using Super-Gaussian pulses or Constant-S* surfaces 

In this subsection we would like to present the results of some numerical experiments. These experiments were 
carried out with the ID version of the VLPL (Virtual Laser Plasma Laboratory) PIC (Particles In Cells) code in 



combination with a 3D cylindrical geometry numerical propagator based on Eq. ( 60 I to obtain the far-field on the 
optical axis from each series of ID simulations. The simulations were made with HHG at a planar surface in mind, 
but in subsection | VI D | we are going to see that all results can easily be re-interpreted to suit HHG at a curved surface 
in confocal geometry. This is crucial in regard of the exciting possibilities opened up by coherent focusing of the 
harmonics radiation. 

We start with a simulation of a Gaussian laser pulse, given by a{x, r, t) = 
Oq Re [exp LLio(a;/c — t) — (t/r)^ — (r/cr)^)] , normally incident onto a surface with a steep and constant dens- 
ity profile. The laser parameters are: ao = 20, a = 5A, r = 27r/a;o. The surface density is = 90 ric, so that the 
ultra-relativistic similarity parameter P6 5*0 — {aQUc) = 4.5 at the maximum of the laser pulse. 



The resulting spectra for the near and the far-field on the optical axis in comparison to the analytical estimate (61 1 



disregarding the phase variation are shown in Fig. 30 a). It is seen, that there is still a big difference in the slope of the 
spectrum compared to Eq. (61 1. Although the spectrum in the far field decays slightly slower than the spectrum close 
to the surface, the improvement is far behind from what we could expect if the attosecond phase remained constant 
along the surface. 



In Fig. 31 we see what happens at an intermediate distance from the surface. According to Eq. (62 1, the pulse 
should be self-focused at a distance of around Xgf sa 130 A. Indeed, we observe that the pulse is self-focused, and 
the self-focusing length is in reasonable agreement with the analytical estimate. The achieved intensity is about four 
times the laser intensity, and it is reached at a distance of about IGOA from the surface. Even attosecond peaks can 
be seen, yet the contrast ratio is quite poor. 

In order to improve the quality of these pulses, we should aim to keep the attosecond phase (p - and thus the 
relativistic ^-parameter - constant alongside the surface. Basically, there are two possibilities to achieve this: 
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Figure 30: Spectra, near- and far-field {x — oo) for different types of laser pulses, normalized to the intensity of the fundamental. 
(ao = 20) 



1. The use of laser pulses with sufficiently flat intensity distributions across the focal spot, e.g. Super-Gaussian. 
This causes the part of the laser pulse that contributes considerably to the HHG to be at a nearly constant 
intensity level. 

2. Varying the surface density in a way, so that S{r) — n^(r) / (a(r) n^) — const. 

For testing these ideas, we perform additional simulation runs: some of them using Super-Gaussian laser pulses 
/ oc exp [— or / oc exp [— (f/cr)*] and some of them using a conventional Gaussian laser pulse, but a surface 
with radially varied density so that the local similarity parameter remains constant S{r) = 4.5. In the simulations 
with the Super-Gaussian laser profile, a was chosen in a way so that the laser power and amplitude are the same as 
in the corresponding simulations with the Gaussian pulse. 



The spectra obtained from these simulations are depicted in Fig. 30 'b)-(d). We see a great improvement compared 



to the unoptimized case Fig. 30 a). Evidently, the spectra decay much slower in t he far- field, reaching close to the 
ideal / cx lo~'^/^ line that was expected from the analytical estimate in subsection VIA The lower frequencies are 
filtered out by diffraction in space. Therefore we may refer to these schemes as "spatial spectral filters". 

The advantage of the constant-S* surface can nicely be seen in Fig. [32] showing the integrand from which the far- 



field is calculated. The integration to obtain the far-field is carried out along the path i! = t — \J 



'-jc. Thus 



for the very far field in front of a planar surface, the integration path becomes a straight line. The intensity of the 
integral becomes maximal if the integrand is big over the whole integration path. Therefore we see, that self-focused 
or defocused attosecond pulses are represented by curved lines in our diagram, whereas non-self-focused attosecond 
pulses are represented by straight lines. Further, longer lines lead to higher peak intensities. With this knowledge, the 
advantage of a constant-^ surface can be easily understood from Fig. |32j Notice also that for the constant-S* surface 
the side peaks are strongly defocused, yielding a better contrast ratio to the main peaks. For Super-Gaussian pulses 
the image would look similar to Fig. [32ja), except for that the upper part is stretched. 

To liven up the picture of how the attosecond pulses emerge, we take a look at the on-axis field at different distances 
X from the surface. For Fig. |33] we choose a simulation with a constant-S* surface, because the process of "vacuum" 
attosecond pulse generation is most pronounced here. While we depart from the surface together with the reflected 
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Figure 31: Self-focusing of a Gaussian pulse (ao — 20, a = 5A) in front of a planar surface compared to focusing of a Super- 
Gaussian pulse (same amplitude and power) to the same distance using a curved surface. The time unit To = 27r/ci;o is just the 
laser period. 



radiation, we see how the attosecond pulses get rectified and the whole rest of the radiation is simply diffracted away 
from the optical axis. 

Focusing these improved pulses using a confocal setting yields a much better result than the self-focusing of a 
Gaussian pulse in front of a planar surface as can be seen from Fig. |3l] Here a Super-Gaussian pulse was chosen, but 
the use of a constant-S* surface leads to a similar effect, as shown further below. 

Let us now have a look at a broader range of parameters. In Fig. [34] the intensity and duration of the attosecond 
pulses in the far-field is compared for different laser amplitudes oq and all the proposed schemes. As shown in the 
following subsection |VID I the results can be applied for planar surfaces as well as for focusing geometries. We compare 
the attosecond pulses in the far-field of a planar surface or at the focal spot in front of a spherically curved surface. We 
notice once again, that the Super-Gaussian laser pulse focal spots and the constant-S" surfaces yield a clear advantage 
for the attosecond pulse generation. 

In earlier works |5J |35], transmission filtering has been suggested as a technique to improve the quality of the 
attosecond pulses. Let us compare this to our method of spatial spectral filtering via shaping of the laser pulse focal 
spot. 

Applying transmission filters directly to the results of ID PIC calculations, as it has been done in previous works, 
yields a somewhat unrealistic picture, since the filters have to be placed inside the far-field in a real experiment. In 
this work we consider the 3D geometry and apply optical filters to the far-field radiation. 

We compare attosecond pulses generated by lasers with a Gaussian and a Super-Gaussian focal spot, see Fig. [35] The 
first thing to notice is that filtering infiuences the temporal structure of the attosecond pulses. While the attosecond 
pulses in the unfiltered far-field are pure half-cycle pulses, optical filtering can generate single- or multi-cycle pulses. 



depending on the filter frequency, see also Sec. IV A Then, unlike transverse pulse shaping (see Fig. 34), transmission 
filtering naturally leads to a decrease of the attosecond peak intensity. To obtain a significantly shortened pulse, one 
needs to use filters with a very high threshold frequency, eating up most of the pulse energy. Nevertheless, filtering 
leads to an improvement of the contrast ratio by a factor of about 3 for totfi = IOOwq- 

To get the shortest possible pulse duration and the best contrast ratio, we recommend to combine transverse pulse 
shaping with the use of an optical filter. When attosecond pulses with a maximum peak intensity are required, the 
transverse pulse shaping or constant-S* surfaces combined with the proper focusing geometry are the best option. 



D. Optical scalings for harmonics focusing 

In this subsection, we examine focusing geometries for the surface harmonics radiation. In Ref. [33] it has been 
shown, that coherent harmonic focusing (CHF) has the potential to produce unprecedentedly intense electromagnetic 
fields. The created intensities may be so extreme, that they can be used to explore exotic QED effects such as vacuum 
polarization or even electron-positron pair creation |43} I44| . However, the conditions under which the phenomenon 
was examined in Ref. [39 were strongly idealized: A perfect spherical wave, uniformly illuminating a curved plasma 
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Figure 32: The integrand ^{r,t) = r dtF-{x — 0,r,t) from Eq. (60 1. The colour scale is the same for both plots. Gaussian 
laser pulses with ao = 20 are used. 



Figure 33: The reflected field on the optical axis, observed at different distances x from the surface. Data from PIC simulation 
using a Gaussian laser pulse with oq = 20 focused on a planar constant-iS surface. 

surface with a large solid angle offl = 1 and a tiny radius of R — AX was studied - something that is not achievable with 
a focused Gaussian beam. Therefore, important effects such as the variation of the laser intensity on the harmonics 
generating surface were not taken into account. Here, we discuss CHF under more realistic conditions. 

In order to better understand CHF, we start by assembling some optical scaling laws for the broadband harmonics 
radiation. These laws are immediate consequences of Eq. (59 1. As the geometries and frequency spectra involved in 



CHF may be unusual, we also make an effort to give conditions of validity for the scaling laws in the cases when they 
are different from the ones for the fundamental Kirchhoff integral ( 59 1 . 



Before we start considering curved surfaces, we have a look at what happens, when the size of the focal spot is 
changed on a planar surface. First, we consider the case when the laser field is focused onto a planar surface and the 
focal spot size is varied but the maximum amplitude of the vector potential ao is kept constant: 

f!')(0, r, t) =f1"^(0, ar, i), 

wherein a is the dimensionless factor describing the focal spot scaling. Now assuming na^ / {4Xx^) <C 1, we can 
calculate that the reflected radiation scales like: 



f'l\x, 0, t-x/c) 



^ F^_^\a'x, Q,t-a^x/c). 



(66) 



Therefore, as long as Eq. (64 1 holds, a variation in the focal spot size will yield an exactly similar field structure in 
the far-field, just scaled in size. 
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Figure 34: Comparison of the intensity and duration (full width half maximum) of the far-field attosecond pulses for different 
laser intensities, radial laser shapes and plasma surfaces. The intensity unit Jo — {aoxm/x)^ is the peak intensity value we 
would get at the same distance from the surface for a Gaussian laser pulse with the same power and maximum vector amplitude 
without the generation of harmonics. The time unit To = 2-k/ujo is the laser period. The constan t g- line (turquoise) in panel 
(a) also denotes the intensity gain that can be expected from the focusing scheme depicted in Fig. 37 Here, 5* = 4.5 was used. 



Now we are interested in actively focusing the harmonics radiation. The most straightforward way to do so is the 
use of a curved surface for HHG, since it can do without relying on any optical components suitable for the extremely 
broadband radiation. We consider a confocal geometry, in the sense that the field distribution F_ of the radiation 

on the curved surface is the same as for the focal spot in the planar geometry p'f"'^ - in both phase and intensity. 
This can be written as: 



F^"\x = - r2 - R, r, t) = f1°°^(0, r, t), 

where R is the radius of the curved surface. We find that the field in front of the planar surface at the distance x, 
is similar to that in front of the curved surface at the distance xn, where Xco and xr are related according to 



Xr Xoo R 



In this case, the field relation is 



F^''\xR,0,t-XR/c) = ^F^^\x^,0,t-x^/c). 

Xr 

This result becomes exact for small focusing solid angles, but also in the most interesting limit xr — > R, which 
corresponds to the actual focal spot when there is no self-focusing present. The scaling law allows us to simply 
re-interpret all results obtained for a planar surface in subsection |VI C| to such for a spherical; of course, as long as 



the condition (64) holds. 

Next, we consider the variation of the focal distance, but keep the field amplitude at the surface and the solid angle 
constant: 

F'f^ar, 0, t) = Fl^^(r, 9, t) 
^ f[°'\0, t - ar/c) = ai^l^'(0, t ~ r/c). (69) 
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Figure 35: The temporal structure of the attosecond pulses in the far-field of a Supergaussian (left column) and a Gaussian 
(right column) laser-pulse, both with the same power and peak amplitude ao = 50 in the focus. The intensity is normalized so 
that the maximum in case of no harmonics and no absorption for the Gaussian Laser would be at I /lo = 1. In the second and 
third lines high-pass filters are applied. The filter function is a simple step like function with a linear transition of the width 
6 = tJo, the threshold frequency is uuth- 



So if the focal distance is varied, the intensity in the focal spot increases proportionally to the input power. Since a 
higher intensity at the surface generally creates a bigger number of harmonics and therefore leads to more than linear 
amplification of the radiation in the focal spot, we should aim for a focal distance as small as possible. 

If the solid angle of the mirroring surface is varied, but the field amplitude and the focal distance are kept constant, 
the intensity in the focal spot increases stronger than the input power: 



Fi"^(i?, a0, t) = f':^\r, 6, t) 
^f1"^(0, t) = a^F^^\0,t). (70) 

This intensity gain is even stronger than the one achieved by the temporal focusing gain due to HHG. Therefore, to 
maximize the focal spot intensity with a constant laser power, it is first needed to maximize the solid angle of CHF. 
Of course, this solid angle is limited by the focusing geometry of the driving laser itself. Second, the CHF distance 
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Figure 36: Illustration of the geometry underlying Eqs. (671, (68 1 
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Figure 37: Suggested geometry for coherent focusing of the surface harmonics radiation, produced by a Gaussian laser beam. 
Ideally, the plasma surface should be designed in a way, that the S-parameter does not change with the distance from the 
optical axis. 



should be minimized in order to make maximum use of the temporal focusing gain due to surface HHG. As shown in 
subsection |VI C[ the use of a constant-5 surface would be ideal here, compare Fig. [34] 

A focusing geometry that fulfils the above criteria is presented in Fig. [37] The laser is focused to a distance of one 
Rayleigh length xr hi front of the HHG surface. The plasma surface is spherically curved with a radius of 2xr and 
the density of the surface is modulated in a way, that the S'-parameter S — n,,/ [aQjic) is constant everywhere. The 
harmonics radiation will then be coherently focused to a distance of 2xfl in front of the curved surface. 

The intensity gain achieved by this scheme is directly related to the idealized gain displayed in Fig. |34j However 
due to the inevitable finite distances from the focus it is reduced to 



1 



ao 



(71) 



Table |III| gives an overview of the gain to be expected from the scheme for different intensities and surface S 
parameters. Note that in the case of S' = 3 and oq = 1000, as can be expected for the European ELI project |46j, an 



(Wcm ^ /im^) in focus 


ao on surface 


'^CHF 






S = 2 


5 = 3 


5 = 4 


2.74 X 10^° 


10 


4.4 


5.3 


3.3 


2.74 X 10^^ 


100 


337 


346 


187 


2.74 X 10^* 


1000 


1185 


5250 


5311 



Table III: Expected gain in temporal peak intensity for coherent focusing of harmonics according to the scheme Fig. |37[ as 
computed by series of ID-PIC simulations in combination with a diffraction integral. 
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output intensity of /(^jjp ~ lO^^Wcm"^ is computed. This would be sufficient to access an entirely new regime of 
physics in which the vacuum itself becomes non-linear due to exotic QED effects. 

In this section, we have discussed the effect of diffraction on the HHG radiation under realistic conditions. To 
make maximum use of these, techniques to shape the focal spot or the plasma surface have to be implemented. 
With these techniques, it may be possible to produce intensities that are more than a thousand times higher than 
the conventional focusing intensity of the laser, opening up the possibility of verifying vacuum QED effects with 
ultraintense laser systems. 



VII. CONCLUSIONS 



The reflection of relativistic light at overdense plasma surfaces, performing a strongly non- linear oscillation, is 
currently one of the most promising candidates for the production of intense attosecond pulses. 

We have taken a fresh look at the theory of their generation. The foundations of the supposedly well-known ROM 
model have been re-investigated, yielding a clearer picture of the scope of application of the model. Further, the 
model has been extended to higher order 7-spikes, demonstrating the possibility of modulated spectral structures 
and power law spectra / cx w^'' with exponents g < 8/3 even within the ROM model. The explicit formulation of 
the TROM model and its comparison to the ROM model sheds additional light on the physics of relativistic high 
harmonics generation. 

A third model was motivated by numerical observations: amazingly dense and narrow electron "nanobunches" 
may form at the plasma surface, emitting coherent synchrotron radiation efficiently. This nanobunching regime of 
relativistic HHG is optimal for attosecond pulse generation in the sense that the generated pulses bear almost the full 
energy of the entire optical cycle of the driving laser. Here, we expect a flat power law spectrum with q < 4/3 up 
to a smooth cut-off at a frequency which is determined either by the bunch relativistic energy a;,.^ cx loqI^ or by the 
nanobunch width uJrf on c/S . 

We have also studied carefully the phase properties of the relativistic harmonics, something that has largely been 
neglected so far. The relation of the electron surface motion to the spectral line structure has been investigated. It 
has been found that the spectral line structure can deliver valuable information on the motion of the electron surface 
on a femtosecond timescale. 

The phase properties are also crucial when considering the free space propagation of the harmonic radiation. We 
have shown that diffraction can be harnessed as a spatial spectral filter for the harmonics radiation by designing the 
target surface or the laser pulse focal spot in the right way. This way, attosecond pulses can be extracted efficiently 
without the use of optical transmission filters. 

Our study can further provide the basis for focusing schemes of the harmonics radiation. CHF (coherent harmonic 
focusing) has the potential to produce intensities that exceed the one of the driving laser by several orders of magnitude. 
Our proposed focusing scheme was shown to yield an intensity enhancement by a factor of 3000 with parameters as 
expected for the ELI facility. Here, unlike previous studies, we took a realistic Gaussian laser pulse geometry into 
account. 



Appendix A: Stationary Phase Method 

In this section, the asymptotic evaluation of integrals via the stationary phase method is explained. The method 
has e.g. been applied to diffraction integrals and the calculation of synchrotron spectra and plays a vital role in the 
theory of high harmonics generation at overdense plasma surfaces. For a comprehensive introduction to this method 
and related ones, consider e.g. the book by Wong |47|. Here, we describe the method briefly with the applications 
from section Hill in mind. 

We are interested in integrals of the form 

/oo 
g{t)exp[iujf{t)] dt, (Al) 
-00 

where g{t) and f{t) are assumed to be smooth functions. We want to find an asymptotic approximation for F{lu) in 
the limit of big uj. Then, the rapidly oscillating integrand cancels everywhere except for the regions of "stationary 
phase" where df{t)/dt « 0. 

In the simplest case, we can find a set of well separated points {tk} on the real axis where df{tk)/dt = and 
d^ f {tk) / dt^ 7^ 0- Then f{t) and g{t) can be Taylor expanded around these points: f{t w tk) ~ ak + bk{t — tkY and 
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g(t « tk) w Cfe. Now. the integral can be evaluated analytically: 

/oo 
exp[iiubk{t-tkf] dt 

k 

It can be seen that the behaviour of the integral depends on two factors. The first factor scales cx oj"^/^ and 
constitutes a spectral envelope, whereas the sum determines if the contributions from each stationary point interfere 
positively or negatively and is thus responsible for the structure of the single "harmonic" lines. This sort of behaviour 
also extends to the more complicated cases discussed below, but the envelope factor varies sensitively according to 
the exact structure of the stationary phase points. Let us now go on to discuss the relevant cases. 



1. First order 7-spikes 



Taking a look at the integrals that we encounter in the models presented in Section |III| we find that none of them 
ever contains points where the condition df{t)/dt = is exactly fulfilled. It is possible to understand this in terms of 
physics. The phase functions f[t) in the models are always connected to the difference between the trajectory of a 
point, that is somehow connected to the plasma motion and the motion of the emitted light wave. Since the plasma 
cannot be faster than light, f{t) is strictly monotonic, consequently df /dt 7^ 0. 

Still it is possible to apply the stationary phase method by considering points where df{t)/dt w 0. Technically, this 
can be viewed as a region where two saddle points, that are located in the complex plane slightly off the real axis, 
closely merge. 

The ROM model, as discussed in Sec. |IIIB] leads us to integrals of the sort 



^±(^) I ey.^[i{t{uj±ujQ) + x{t){uj^uja))] {I + x{t)) dt. 



(A2) 

This integral has to be handled with attention: Note that although w ^ is assumed, we must not neglect in 
the exponent. The reason for this will become evident later. 

The Taylor expansion of x{t) around the velocity maximum can be written down as x (t) ^ -vt + at^/3. In the 
case a 7^ 0, we speak about a "7-spike of the order 1". In the case, when a — 0, higher orders of the Taylor expansion 
have to be considered. It is discussed in subsection lA 21 

We shift the stationary phase point to to = without loss of generality here, as we are not interested in absolute 
phase terms. Using the abbreviation S = 1 — v ^ 1/(27^), we get F{uj) = Fi{uj) + F2{uj) with 



Fi{uj) = S J exp 
F2{io) = a I t^ cxp 



i t {5u) ± (2 — 5)ujo) + a{uj ^ Wq) 



dt 



I t {5uj ± (2 — + aiu ^ wo)" 



(A3) 



dt 



By now it should become clear why F2 and the above mentioned wo-terms could not be neglected: We want to presume 
Lo 3> wqi but not 5uj 3> wq- Later on, we may neglect loq compared to u and 5 compared to 1. 

-Fi and F2 can now be expressed in terms of the well-known Airy function Ai(a;) = (27r)~^ ^-'^P (* {^^ I'^T) 



Fdu;) = f^'^ AiiO (A4) 
{/a (w =F Wo) 

F2{u;) = ^^Ai(0. (A5) 
w =F Wo 

with ^ = {Slu ± (2 — S)ujo) / (a (w ^ ujq))^^^. For the calculation of F2, we made use of Ai"(a;) = a;Ai(a;). After taking 
the sum of Fi and F2, the w^^/"^ terms cancel and only the uj^^^^ term remains, which represents the leading order 



now: 



F±H^ .4/3 Ai(0- (A6) 



3/ 



a (w =F Wo 
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Taking the absolute square yields the famous — 8/3-power law spectrum. 
For the TROM model (section IIIC), the integral looks a bit different: 



(A7) 

The integration works in complete analogy to the case shown above, but in this case, the ^"^/^-ternis do not cancel 
out. Therefore, we obtain in highest order: 



An 



{/a (w=Fwo) 



Ai(0 



(A8) 



In the case of CSE (section HID I, the integral is of the sort: 



F{uj) = J y{t) exp [-iuj (t + x{t))] dt. (A9) 

To get some meaningful result out of this, we need to make an assumption about the relation between y{t) and 
x{t). We assume, that during the time of harmonic generation, the absolute velocity {x^ + y^Y^"^ is approximately 
constant and close to the speed of light. This is reasonable in the ultra-relativistic regime. With this assumption, the 
stationary phase points are exactly the points, where y vanishes and the electrons move towards the observer. Now 
we can Taylor expand y{t) ~ a^t and x{t) ~ —vt + ait^ /i. Substituting into Eq. (A9l yields: 



dt. 



F{lu) — ao J t exp i ^—coSt — ioai — 
where S ^ 1 — v as above. Again, the result can be expressed in terms of the Airy function: 

F(^) 



-2TTaQi , . / Jo;-^/-^ 
Ai' ' 



(aiw) 



2/3 



(AlO) 



(All) 



2. Higher order 7-spikes 



In the previous subsection (Sec. All, we dealt with the case when the transverse velocity of the electrons pass zero. 
Together with the assumption of ultrarelativistic motion this lead us to the Taylor expansion x{t) = —vt + at^ /3. In 
this section we deal with the possibility, that the transverse velocity does not go through, but touches zero, so that 
the third order of x{t) vanishes. 

In general, if the first 2n orders of x{t) vanish, it can be written: x{t) — —vt + at^"+^/(2n + 1). We refer to this 
case as a "7-spike of the order n". Inserting this into Eq. (A2| yields: 



F2{0J) 



exp 



i t {5uj ± (2 — S)ujo) + a{uj ^ ujq) 



2n+l 



a I t^" exp 



2n + l 

i\t{5Lo±{2- 5)L0f)) + a{Lo =F ^o) — — r 

2n + 1 



dt 



(A12) 



dt 



These integrals can now be expressed by a generalized Airy function, which we define as gAi„(a:;) — 
(27r)~^/^exp \i (xi + i^"+^/(2n + 1))] dt. Note that for n = 1 we retain the Airy function and for n = 2 we 
obtain a special case of the canonical swallowtail integral j29j . Since the gAi„ (a;) are not available in general purpose 
numerical libraries, their numerical evaluation is explained in Sec. [C| 

In analogy to the Airy function, the gAi„(x) fulfil ODEs: d^"gAi„(x)/c?a::^" + (— l)"a; gAi„(a::) — 0. Exploiting this, 
Fi and F2 become 



Fi(w) 
F2(w) 



2ti5 



(A13) 
(A14) 
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where f = {Slo ± (2 — 6)u}o) / (a (w =F wq))^^*'^"^^''. After again taking the sum of Fi and F2, the w terms 
cancel and what remains is: 



±47r 



(A15) 



In complete analogy the TROM model, represented by the integral (A7l, yields: 



47r 



?Ai„(e)- 



(A16) 



Finally we calculate the CSE integral (A9l for arbitrary orders of the 7-spike. Here, this means y{t) — a^V^ and, 
consequently, x{t) = —vt + Qi<^"+^/(2n + 1). Since this works in complete analogy to the hitherto discussed cases, 
we just present the result: 



-27^«o^" rf"gAi„(e) 
(aia;)"+V2"+i d^n ' 



(A17) 



wherein ^ = c^2"/2n+i^/Q,i/(2»+i). 



Appendix B: Lorentz Transformation to Describe Oblique Incidence in ID 

In laser-plasma theory, it is often convenient to describe things in a one dimensional (ID) slab geometry, i.e. all 
spatial derivatives perpendicular to the optical axis are neglected. For not too small laser focal spot sizes, this is 
very often a reasonable approximation and leads to great simplifications in numerical as well as analytical theory. In 
numerical computations, the grid size can be reduced by orders of magnitude, allowing for higher resolution in the 
critical dimension. In analytical calculations, it sometimes enables us to give closed form solutions and straightforward, 
comprehensible models. Whereas it is obvious that the ID treatment can be employed in situations of normal laser 
incidence, it can also be extended to oblique laser incidence. Therefore, as shown by Bourdier in Ref. [48 , a Lorentz 
transformation does the job. 

From the lab frame £, we transform to the inertial frame S, in which the laser is normally incident. Let the 
light wave in C be described by the frequency 4-vector (cjg , ck^ , cky, O) = lJq (1, cos a, sin a, 0), wherein x denotes 

the direction normal to the surface and a is the angle of incidence. In 5, we claim cky = 0. Thus, the Lorentz 
transformation is given by the matrix 



A = 



/ 7 -/?7 0\ 

10 

-/37 7 

V 1 / 



(Bl) 



with l3 — sin a and 7 = (cos a) ^. Now, we can derive all interesting magnitudes. In the frame S, the plasma is 
streaming with the velocity Vy — —cf3, the laser wavelength is altered by A"^ = jX^ and therefore the corresponding 
critical density changes to nfj. = 7^^nf^. The electron density itself changes to = in^ , thus the normalized density 



n" 



(B2) 



scales even with 7"^ 



Appendix C: Numerical Calculation of the Generalized Airy Function 

This section explains the numerical computation of the integral 



1 



.Ai„(x)^- / e\^'^^)dt. 
It: 



(CI) 
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The unmindful immediate application of a trapezoidal formula would fail here due to the rapidly oscillating in- 
tegrand, which further does not vanish at infinity. The integral only converges because of the steadily decreasing 
oscillation period for t — ?> ±00. 

However, a simple trick can be applied to calculate the integral numerically: We shift the contour of integration 
along the imaginary axis by a margin of a > 0. Due to Cauchy's integral theorem, this will not change the results, 
since for large t the contours connecting the real axis to the new integration path do not contribute. Next, the 
symmetry of the integrand can be exploited, so that the contour of integration can be halved. We obtain: 

gAiJx) = - / cxp[-Im(/(x,i))] xcos[Re(/(x,t))] (C2) 
^ Jo 

(t + ml2"+i 

f{x,t) = + ^ ; ^ ■ (C3) 

2n + 1 

Note that, since a > 0, f{x, t) possesses the highly desirable property limt^oo [^"^(fi^it))] = +oo. This means. 



that the integrand vanishes exponentially for large t, and we can approximate (C2) with a trapezoidal formula. 

Since a does not affect the value of the integral, it can be chosen in a way so as to minimize the computational 
effort. For too small a, the integrand vanishes only slowly for t — >■ 00 so that the numerical upper boundary would 
have to be very high. For too large a, the integrand oscillates more rapidly, so that the time step would have to be 
very small. Our experience showed, that the best choice for a depends mainly on the order n of the function. For 
71 < 2 we found a = 1 a good choice, whereas for n > 2, smaller a work better. 
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